compiler: prevent undefined Temp through global init

Author: mloubout

devito/passes/clusters/aliases.py

@@ -12,6 +12,7 @@
     Interval, IntervalGroup, IterationSpace, LabeledVector, Queue, Vector, extrema,
     maximum, minimum, normalize_properties, relax_properties, unbounded, vmax, vmin
 )
+from devito.ir.support import null_ispace
 from devito.passes.clusters.cse import _cse
 from devito.symbolics import (
     Uxmapper, estimate_cost, retrieve_functions, reuse_if_untouched, search, sympy_dtype,
@@ -131,7 +132,8 @@ def _aliases_from_clusters(self, cgroup, exclude, meta):
             # AliasList -> Schedule
             schedule = lower_aliases(aliases, meta, self.opt_maxpar)
 
-            variants.append(Variant(schedule, exprs))
+            if schedule:
+                variants.append(Variant(schedule, exprs))
 
         if not variants:
             return []
@@ -860,6 +862,7 @@ def lower_schedule(schedule, meta, sregistry, opt_ftemps, opt_min_dtype,
     make = TempFunction if opt_ftemps else TempArray
 
     clusters = []
+    inits = []
     subs = {}
     for pivot, writeto, ispace, aliaseds, indicess in schedule:
         name = sregistry.make_name()
@@ -928,8 +931,11 @@ def lower_schedule(schedule, meta, sregistry, opt_ftemps, opt_min_dtype,
             assert writeto.size == 0
 
             dtype = sympy_dtype(pivot, base=meta.dtype, smin=opt_min_dtype)
-            obj = Temp(name=name, dtype=dtype)
+            const = not meta.guards
+            obj = Temp(name=name, dtype=dtype, is_const=const)
             expression = Eq(obj, uxreplace(pivot, subs))
+            if not const:
+                inits.append(Eq(obj, 0))
 
             callback = lambda idx: obj  # noqa: B023
 
@@ -959,6 +965,13 @@ def lower_schedule(schedule, meta, sregistry, opt_ftemps, opt_min_dtype,
         # Finally, build the alias Cluster
         clusters.append(Cluster(expression, ispace, meta.guards, properties))
 
+    if inits:
+        # To avoid undefined variables when an constant (Temp) alias is used
+        # within different guards/loop, we need to initialize it outside of the loops
+        # so that it's globally defined.
+        # See tests/test_dse.py::TestAliases::test_split_cond
+        clusters.insert(0, Cluster(inits, null_ispace))
+
     return clusters, subs
 
 

devito/passes/clusters/cse.py

@@ -400,6 +400,10 @@ def _(exprs):
 
 @_catch.register(sympy.Eq)
 def _(expr):
+    if isinstance(expr.lhs, Temp):
+        # Prevent re-factorization of constant alias
+        return {}
+
     mapper = _catch(expr.rhs)
     try:
         cond = expr.conditionals

devito/types/basic.py

@@ -1589,7 +1589,7 @@ def _rebuild(self, *args, **kwargs):
         comps = [f.func(*args, name=f.name.replace(self.name, newname), **kwargs)
                  for f in self.flat()]
         # Rebuild the matrix with the new components
-        return self._new(comps)
+        return self._new(*self.shape, comps)
 
     func = _rebuild
 

examples/mpi/overview.ipynb

@@ -460,7 +460,7 @@
        "  _MM_SET_DENORMALS_ZERO_MODE(_MM_DENORMALS_ZERO_ON);\n",
        "  _MM_SET_FLUSH_ZERO_MODE(_MM_FLUSH_ZERO_ON);\n",
        "\n",
-       "  float r0 = 1.0F/h_x;\n",
+       "  const float r0 = 1.0F/h_x;\n",
        "\n",
        "  for (int time = time_m, t0 = (time)%(2), t1 = (time + 1)%(2); time <= time_M; time += 1, t0 = (time)%(2), t1 = (time + 1)%(2))\n",
        "  {\n",

tests/test_dse.py

@@ -2685,6 +2685,25 @@ def test_space_and_time_invariant_together(self):
             'tx0_blk0y0_blk0xyzyz'
         )
 
+    def test_split_cond(self):
+        grid = Grid((11, 11))
+        time = grid.time_dim
+
+        u = TimeFunction(name='u', grid=grid, time_order=2, space_order=2)
+
+        ct = ConditionalDimension(name='ct', parent=time, factor=2)
+        ct2 = ConditionalDimension(name='ct2', parent=time, factor=4)
+
+        eq0 = Eq(u.forward, u + cos(time), implicit_dims=ct)
+        eq1 = Eq(u.forward, u.forward + 1, implicit_dims=ct2)
+        eq2 = Eq(u.forward, u.forward + cos(time), implicit_dims=ct)
+
+        op = Operator([eq0, eq1, eq2])
+        cond = FindNodes(Conditional).visit(op)
+        assert len(cond) == 3
+        assert str(cond[0].args['then_body'][0].exprs[0]) == 'r0 = cos(time);'
+        assert str(op.body.body[0].body[0].body[0]) == 'float r0 = 0;'
+
 
 class TestIsoAcoustic: