compiler: Avoid CSE capture of opaque DefFunction calls

devito/passes/clusters/cse.py

@@ -13,7 +13,9 @@
 
 from devito.finite_differences.differentiable import IndexDerivative
 from devito.ir import Cluster, Scope, cluster_pass
-from devito.symbolics import Reserved, estimate_cost, q_leaf, q_terminal, search
+from devito.symbolics import (
+    DefFunction, Reserved, estimate_cost, q_leaf, q_terminal, search
+)
 from devito.symbolics.manipulation import _uxreplace
 from devito.tools import DAG, as_list, as_tuple, extract_dtype, frozendict
 from devito.types import Eq, Symbol, Temp
@@ -462,3 +464,14 @@ def _(expr):
     mapper[Candidate(expr)].append(expr)
 
     return mapper
+
+
+@_catch.register(DefFunction)
+def _(expr):
+    """
+    Handler for opaque C-level calls (e.g. `make_float4(...)`). Their return
+    dtype is invisible to `extract_dtype` -- capturing one would bind it to a
+    temporary of the wrong type (`float r0 = make_float4(...)`) -- and they
+    are not guaranteed to be pure, so they are left untouched.
+    """
+    return {}

tests/test_cse.py

@@ -11,7 +11,7 @@
 from devito.ir import Conditional, DummyEq, FindNodes, FindSymbols
 from devito.ir.support import generator
 from devito.passes.clusters.cse import CTemp, _cse
-from devito.symbolics import indexify
+from devito.symbolics import DefFunction, indexify
 from devito.types import Array, Symbol, Temp
 
 
@@ -131,6 +131,29 @@ def test_temp_order():
     assert type(args[2]) is CTemp
 
 
+def test_deffunction_not_captured():
+    """
+    Opaque C-level calls (DefFunction) must not be CSE-captured: their return
+    dtype is invisible to the type inference (e.g. `make_float4` returns a
+    `float4`, but `extract_dtype` only sees the `float` arguments), so binding
+    one to a temporary declares it with the wrong type and the generated code
+    does not compile.
+    """
+    t0 = Symbol(name='t0', dtype=np.float32)
+    t1 = Symbol(name='t1', dtype=np.float32)
+    a = CTemp(name='a', dtype=np.float32)
+    b = CTemp(name='b', dtype=np.float32)
+    call = DefFunction('make_float4', (t0 + t1, t0 + t1, t0, t1))
+    exprs = [DummyEq(a, call), DummyEq(b, call)]
+
+    counter = generator()
+    make = lambda _: CTemp(name=f'r{counter()}')
+    processed = _cse(exprs, make)
+
+    # The two calls are left in place, uncaptured
+    assert processed == exprs
+
+
 def test_w_conditionals():
     grid = Grid(shape=(10, 10, 10))
     x, _, _ = grid.dimensions