Andrew add tests to test operator 2194

.github/workflows/asv.yml

@@ -58,7 +58,7 @@ jobs:
 
     - name: Run benchmarks
       run: |
-        asv run -v --strict --show-stderr --config benchmarks/regression/asv.conf.json --cpu-affinity 0-7 --machine i7-6700K
+        asv run -v --show-stderr --config benchmarks/regression/asv.conf.json --cpu-affinity 0-7 --machine i7-6700K
 
     - name: Checkout asv-results branch
       uses: actions/checkout@v3

devito/ir/support/basic.py

@@ -426,24 +426,31 @@ def distance(self, other):
             else:
                 v = i - j
                 if v.is_Number and v.is_finite:
-                    # If i and j are numbers, we append an (I) distance
                     if i.is_Number and j.is_Number:
                         return Vector(S.ImaginaryUnit)
-                    # If both i and j are not numbers, there may be dimension-dependent
-                    # dependencies so we append the distance
                     else:
+                        # If both i and j are not numbers,
+                        # there may be dimension-dependent dependencies
+                        # so we append the distance
+
+                        # For example:
+                        # Eq(u[0, y], 1)
+                        # Eq(u[1, 1], u[0, y+1])
                         ret.append(v)
 
-                # We are writing over an entire dimension but reading from one point.
-                # If there are overlaps between the two then we would have a dependency
-                # This is a conservative estimation as there are cases (example below)
-                # where we may or may not have a dependency given that we don't write
+                # We are writing over an entire dimension
+                # but reading from one point.
+                # If there are overlaps between the two
+                # then we would have a dependency
+                # This is a conservative estimation as there are cases
+                # where we may or may not have a dependency
+                # given that we don't write
                 # depending on domain size, which is not compilation-time known
 
                 # For example:
                 # Eq(u[0, y], 1)
                 # Eq(u[1, y+1], u[0, 1])
-                elif i.is_Number and not j.is_Number:
+                elif (not i.is_Number or not j.is_Number):
                     ret.append(S.Infinity)
 
         return Vector(*ret)

devito/passes/clusters/factorization.py

@@ -142,13 +142,16 @@ def run(expr):
                     terms.append(i)
 
             # Collect common funcs
-            w_funcs = Add(*w_funcs, evaluate=False)
-            w_funcs = collect(w_funcs, funcs, evaluate=False)
-            try:
-                terms.extend([Mul(k, collect_const(v), evaluate=False)
-                              for k, v in w_funcs.items()])
-            except AttributeError:
-                assert w_funcs == 0
+            if len(w_funcs) > 1:
+                w_funcs = Add(*w_funcs, evaluate=False)
+                w_funcs = collect(w_funcs, funcs, evaluate=False)
+                try:
+                    terms.extend([Mul(k, collect_const(v), evaluate=False)
+                                  for k, v in w_funcs.items()])
+                except AttributeError:
+                    assert w_funcs == 0
+            else:
+                terms.extend(w_funcs)
 
             # Collect common pows
             w_pows = Add(*w_pows, evaluate=False)

devito/symbolics/queries.py

@@ -20,7 +20,8 @@
 # * Number
 # * Symbol
 # * Indexed
-extra_leaves = (FieldFromPointer, FieldFromComposite, IndexedBase, AbstractObject)
+extra_leaves = (FieldFromPointer, FieldFromComposite, IndexedBase, AbstractObject,
+                IndexedPointer)
 
 
 def q_symbol(expr):
@@ -31,7 +32,9 @@ def q_symbol(expr):
 
 
 def q_leaf(expr):
-    return expr.is_Atom or expr.is_Indexed or isinstance(expr, extra_leaves)
+    return (expr.is_Atom or
+            expr.is_Indexed or
+            isinstance(expr, extra_leaves))
 
 
 def q_indexed(expr):
@@ -51,7 +54,7 @@ def q_derivative(expr):
 def q_terminal(expr):
     return (expr.is_Symbol or
             expr.is_Indexed or
-            isinstance(expr, extra_leaves + (IndexedPointer,)))
+            isinstance(expr, extra_leaves))
 
 
 def q_routine(expr):

devito/types/dense.py

@@ -1040,6 +1040,7 @@ def __indices_setup__(cls, *args, **kwargs):
             dimensions = grid.dimensions
 
         if args:
+            assert len(args) == len(dimensions)
             return tuple(dimensions), tuple(args)
 
         # Staggered indices
@@ -1449,16 +1450,10 @@ class SubFunction(Function):
     """
     A Function bound to a "parent" DiscreteFunction.
 
-    A SubFunction hands control of argument binding and halo exchange to its
-    parent DiscreteFunction.
+    A SubFunction hands control of argument binding and halo exchange to the
+    DiscreteFunction it's bound to.
     """
 
-    __rkwargs__ = Function.__rkwargs__ + ('parent',)
-
-    def __init_finalize__(self, *args, **kwargs):
-        super(SubFunction, self).__init_finalize__(*args, **kwargs)
-        self._parent = kwargs['parent']
-
     def __padding_setup__(self, **kwargs):
         # SubFunctions aren't expected to be used in time-consuming loops
         return tuple((0, 0) for i in range(self.ndim))
@@ -1470,12 +1465,8 @@ def _arg_values(self, **kwargs):
         if self.name in kwargs:
             raise RuntimeError("`%s` is a SubFunction, so it can't be assigned "
                                "a value dynamically" % self.name)
-        else:
-            return self._parent._arg_defaults(alias=self._parent).reduce_all()
 
-    @property
-    def parent(self):
-        return self._parent
+        return self._arg_defaults(alias=self)
 
     @property
     def origin(self):