diff --git a/conftest.py b/conftest.py index 65bfddb71c..edf83799ef 100644 --- a/conftest.py +++ b/conftest.py @@ -479,3 +479,9 @@ def check_array(array, exp_halo, exp_shape, rotate=False): assert tuple(array.halo) == exp_halo assert tuple(shape) == tuple(exp_shape) + + +# Main body in Operator IET, depending on ISA +def body0(op): + bidx = 0 if 'sse' not in configuration['platform'].known_isas else 1 + return op.body.body[bidx] diff --git a/devito/passes/clusters/aliases.py b/devito/passes/clusters/aliases.py index a126f96684..5a65b9ce6d 100644 --- a/devito/passes/clusters/aliases.py +++ b/devito/passes/clusters/aliases.py @@ -126,6 +126,8 @@ def _aliases_from_clusters(self, cgroup, exclude, meta): for mapper in self._generate(cgroup, exclude): # Clusters -> AliasList found = collect(mapper.extracted, meta.ispace, self.opt_minstorage) + if not found: + continue exprs, aliases = self._choose(found, cgroup, mapper) # AliasList -> Schedule @@ -271,7 +273,6 @@ def _do_generate(self, exprs, exclude, cbk_search, cbk_compose=None): free_symbols = i.free_symbols if {a.function for a in free_symbols} & exclude: continue - mapper.add(i, make, terms) return mapper @@ -279,6 +280,8 @@ def _do_generate(self, exprs, exclude, cbk_search, cbk_compose=None): class CireInvariants(CireTransformerLegacy, Queue): + _q_guards_in_key = True + def __init__(self, sregistry, options, platform): super().__init__(sregistry, options, platform) @@ -928,7 +931,7 @@ 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) + obj = Temp(name=name, dtype=dtype, is_const=True) expression = Eq(obj, uxreplace(pivot, subs)) callback = lambda idx: obj # noqa: B023 diff --git a/devito/passes/clusters/cse.py b/devito/passes/clusters/cse.py index 9108dadaba..b6ce7c58cd 100644 --- a/devito/passes/clusters/cse.py +++ b/devito/passes/clusters/cse.py @@ -40,6 +40,12 @@ def cse_dtype(exprdtype, cdtype): Return the dtype of a CSE temporary given the dtype of the expression to be captured and the cluster's dtype. """ + if np.issubdtype(cdtype, np.floating) and np.issubdtype(exprdtype, np.integer): + # Integer expression and floating-point cluster: promote to the floating point + # np.promote_types upcast integers (e.g int32 -> Float64) so we + # need to ensure that the promoted type is not larger than the cluster's dtype + return cdtype + if np.issubdtype(cdtype, np.complexfloating): return np.promote_types(exprdtype, cdtype(0).real.__class__).type else: @@ -97,8 +103,9 @@ def cse(cluster, sregistry=None, options=None, **kwargs): if cluster.is_fence: return cluster - make_dtype = lambda e: cse_dtype(e.dtype, dtype) - make = lambda e: CTemp(name=sregistry.make_name(), dtype=make_dtype(e)) + def make(e): + edtype = cse_dtype(e.dtype, dtype) + return CTemp(name=sregistry.make_name(), dtype=edtype) exprs = _cse(cluster, make, min_cost=min_cost, mode=mode) diff --git a/devito/passes/clusters/misc.py b/devito/passes/clusters/misc.py index f248061252..a92c5495f1 100644 --- a/devito/passes/clusters/misc.py +++ b/devito/passes/clusters/misc.py @@ -97,12 +97,12 @@ def callback(self, clusters, prefix): properties = c.properties.filter(key) - # Lifted scalar clusters cannot be guarded - # as they would not be in the scope of the guarded clusters - # unless the guard is for an outer dimension - guards = {} if c.is_scalar and not (prefix[:-1] and c.guards) else c.guards - - lifted.append(c.rebuild(ispace=ispace, properties=properties, guards=guards)) + # If `c` is made of scalar expressions within guards, then we must keep + # it close to the adjacent Clusters for correctness + if c.is_scalar and c.guards and ispace: + processed.append(c.rebuild(ispace=ispace, properties=properties)) + else: + lifted.append(c.rebuild(ispace=ispace, properties=properties)) return lifted + processed diff --git a/devito/symbolics/search.py b/devito/symbolics/search.py index cfb9d504be..f0acdeca72 100644 --- a/devito/symbolics/search.py +++ b/devito/symbolics/search.py @@ -96,7 +96,7 @@ def visit_preorder_first_hit(self, expr: Expression) -> Iterator[Expression]: def search(exprs: Expression | Iterable[Expression], - query: type | Callable[[Any], bool], + query: type | tuple[type, ...] | Callable[[Any], bool], mode: Mode = 'unique', visit: Literal['dfs', 'bfs', 'bfs_first_hit'] = 'dfs', deep: bool = False) -> List | set[Expression]: @@ -104,7 +104,10 @@ def search(exprs: Expression | Iterable[Expression], assert mode in ('all', 'unique'), "Unknown mode" - Q = (lambda obj: isinstance(obj, query)) if isinstance(query, type) else query + if isinstance(query, (type, tuple)): + Q = lambda obj: isinstance(obj, query) + else: + Q = query # Search doesn't actually use a BFS (rather, a preorder DFS), but the terminology # is retained in this function's parameters for backwards compatibility diff --git a/devito/types/basic.py b/devito/types/basic.py index 89de6c5a01..3239638a98 100644 --- a/devito/types/basic.py +++ b/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 diff --git a/examples/mpi/overview.ipynb b/examples/mpi/overview.ipynb index bb7e010779..b46f49aaf4 100644 --- a/examples/mpi/overview.ipynb +++ b/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", diff --git a/examples/performance/00_overview.ipynb b/examples/performance/00_overview.ipynb index 6e80e74d56..72f1836aa4 100644 --- a/examples/performance/00_overview.ipynb +++ b/examples/performance/00_overview.ipynb @@ -508,7 +508,7 @@ "}\n", "STOP(section0,timers)\n", "\n", - "float r1 = 1.0F/h_y;\n", + "const float r1 = 1.0F/h_y;\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", @@ -572,7 +572,13 @@ "+ u[t1][x + 4][y + 4][z + 4] = (f[x + 1][y + 1][z + 1]*f[x + 1][y + 1][z + 1])*((-6.66666667e-1F*r0)*(8.33333333e-2F*r0*u[t0][x + 4][y + 1][z + 4] - 6.66666667e-1F*r0*u[t0][x + 4][y + 2][z + 4] + 6.66666667e-1F*r0*u[t0][x + 4][y + 4][z + 4] - 8.33333333e-2F*r0*u[t0][x + 4][y + 5][z + 4]) + (-8.33333333e-2F*r0)*(8.33333333e-2F*r0*u[t0][x + 4][y + 4][z + 4] - 6.66666667e-1F*r0*u[t0][x + 4][y + 5][z + 4] + 6.66666667e-1F*r0*u[t0][x + 4][y + 7][z + 4] - 8.33333333e-2F*r0*u[t0][x + 4][y + 8][z + 4]) + (8.33333333e-2F*r0)*(8.33333333e-2F*r0*u[t0][x + 4][y][z + 4] - 6.66666667e-1F*r0*u[t0][x + 4][y + 1][z + 4] + 6.66666667e-1F*r0*u[t0][x + 4][y + 3][z + 4] - 8.33333333e-2F*r0*u[t0][x + 4][y + 4][z + 4]) + (6.66666667e-1F*r0)*(8.33333333e-2F*r0*u[t0][x + 4][y + 3][z + 4] - 6.66666667e-1F*r0*u[t0][x + 4][y + 4][z + 4] + 6.66666667e-1F*r0*u[t0][x + 4][y + 6][z + 4] - 8.33333333e-2F*r0*u[t0][x + 4][y + 7][z + 4]))*sinf(f[x + 1][y + 1][z + 1]);\n", " }\n", " }\n", - " }\n", + " }\n" + ] + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ "\n" ] } @@ -652,7 +658,7 @@ }, { "cell_type": "code", - "execution_count": 13, + "execution_count": 12, "metadata": {}, "outputs": [ { @@ -712,7 +718,7 @@ }, { "cell_type": "code", - "execution_count": 14, + "execution_count": 13, "metadata": {}, "outputs": [ { @@ -753,7 +759,14 @@ " }\n", " }\n", " STOP(section0,timers)\n", - "}\n" + "}" + ] + }, + { + "name": "stdout", + "output_type": "stream", + "text": [ + "\n" ] } ], @@ -772,7 +785,7 @@ }, { "cell_type": "code", - "execution_count": 15, + "execution_count": 14, "metadata": {}, "outputs": [ { @@ -863,7 +876,7 @@ }, { "cell_type": "code", - "execution_count": 16, + "execution_count": 15, "metadata": {}, "outputs": [ { @@ -919,7 +932,7 @@ }, { "cell_type": "code", - "execution_count": 17, + "execution_count": 16, "metadata": {}, "outputs": [ { @@ -976,7 +989,7 @@ }, { "cell_type": "code", - "execution_count": 18, + "execution_count": 17, "metadata": {}, "outputs": [], "source": [ @@ -994,7 +1007,7 @@ }, { "cell_type": "code", - "execution_count": 19, + "execution_count": 18, "metadata": {}, "outputs": [ { @@ -1044,7 +1057,7 @@ }, { "cell_type": "code", - "execution_count": 20, + "execution_count": 19, "metadata": {}, "outputs": [ { @@ -1112,7 +1125,7 @@ }, { "cell_type": "code", - "execution_count": 21, + "execution_count": 20, "metadata": {}, "outputs": [ { @@ -1192,7 +1205,7 @@ " }\n", " STOP(section0,timers)\n", "\n", - " float r1 = 1.0F/h_y;\n", + " const float r1 = 1.0F/h_y;\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", @@ -1279,7 +1292,7 @@ }, { "cell_type": "code", - "execution_count": 22, + "execution_count": 21, "metadata": {}, "outputs": [ { @@ -1304,7 +1317,7 @@ "}\n", "STOP(section0,timers)\n", "\n", - "float r1 = 1.0F/h_y;\n", + "const float r1 = 1.0F/h_y;\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", @@ -1369,7 +1382,7 @@ }, { "cell_type": "code", - "execution_count": 23, + "execution_count": 22, "metadata": {}, "outputs": [ { @@ -1404,7 +1417,7 @@ }, { "cell_type": "code", - "execution_count": 24, + "execution_count": 23, "metadata": {}, "outputs": [ { @@ -1483,7 +1496,7 @@ " }\n", " STOP(section0,timers)\n", "\n", - " float r1 = 1.0F/h_y;\n", + " const float r1 = 1.0F/h_y;\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", @@ -1546,7 +1559,7 @@ }, { "cell_type": "code", - "execution_count": 25, + "execution_count": 24, "metadata": {}, "outputs": [ { @@ -1624,8 +1637,8 @@ " }\n", " STOP(section0,timers)\n", "\n", - " float r1 = 1.0F/h_x;\n", - " float r2 = 1.0F/h_y;\n", + " const float r1 = 1.0F/h_x;\n", + " const float r2 = 1.0F/h_y;\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", @@ -1718,7 +1731,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.10.12" + "version": "3.13.11" } }, "nbformat": 4, diff --git a/examples/performance/01_gpu.ipynb b/examples/performance/01_gpu.ipynb index 54eeb2e9fc..3e6e68b57c 100644 --- a/examples/performance/01_gpu.ipynb +++ b/examples/performance/01_gpu.ipynb @@ -142,8 +142,13 @@ "name": "stderr", "output_type": "stream", "text": [ - "NUMA domain count autodetection failed, assuming 1\n", - "Operator `Kernel` ran in 0.01 s\n", + "NUMA domain count autodetection failed, assuming 1\n" + ] + }, + { + "name": "stderr", + "output_type": "stream", + "text": [ "Operator `Kernel` ran in 0.01 s\n" ] } @@ -292,9 +297,9 @@ " const int x_stride0 = x_fsz0*y_fsz0;\n", " const int y_stride0 = y_fsz0;\n", "\n", - " float r0 = 1.0F/dt;\n", - " float r1 = 1.0F/(h_x*h_x);\n", - " float r2 = 1.0F/(h_y*h_y);\n", + " const float r0 = 1.0F/dt;\n", + " const float r1 = 1.0F/(h_x*h_x);\n", + " const float r2 = 1.0F/(h_y*h_y);\n", "\n", " for (int time = time_m; time <= time_M; time += 1)\n", " {\n", @@ -340,7 +345,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.13.5" + "version": "3.13.11" } }, "nbformat": 4, diff --git a/tests/test_dse.py b/tests/test_dse.py index 24eb0cda2b..8399261563 100644 --- a/tests/test_dse.py +++ b/tests/test_dse.py @@ -5,8 +5,8 @@ from sympy import Mul # noqa from conftest import ( # noqa - _R, EVAL, assert_blocking, assert_structure, check_array, get_arrays, get_params, - skipif + _R, EVAL, assert_blocking, assert_structure, body0, check_array, get_arrays, + get_params, skipif ) from devito import ( # noqa NODE, Abs, ConditionalDimension, Constant, DefaultDimension, Derivative, Dimension, @@ -28,6 +28,7 @@ ) from devito.tools import as_tuple from devito.types import PrecomputedSparseTimeFunction, Scalar, Symbol +from devito.types.misc import Temp from examples.seismic import AcquisitionGeometry, demo_model from examples.seismic.acoustic import AcousticWaveSolver from examples.seismic.tti import AnisotropicWaveSolver @@ -347,8 +348,8 @@ def test_scalar_cond(self): trees = retrieve_iteration_tree(op) assert len(trees) == 3 - assert_structure(op, ['t', 't,x,y', 't,x,y'], 'txyxy') - assert trees[0].dimensions == [time] + assert_structure(op, ['t,x,y', 't', 't,x,y'], 'txyxy') + assert trees[1].dimensions == [time] class TestAliases: @@ -2551,6 +2552,7 @@ def test_invariants_with_conditional(self): eqn = Eq(u, u - (cos(time_sub * factor * f) * uf)) op = Operator(eqn, opt='advanced') + assert_structure(op, ['t', 't,fd', 't,fd,x,y'], 't,fd,x,y') # Make sure it compiles _ = op.cfunction @@ -2685,6 +2687,190 @@ 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]) + op(time=5) + + cond = FindNodes(Conditional).visit(op) + assert len(cond) == 3 + # Each guard should have its own alias for cos(time) + assert 'float r0 = cos(time);' in str(body0(op)) + scalars = [i for i in FindSymbols().visit(op) if isinstance(i, Temp)] + assert len(scalars) == 2 + + def test_split_cond_multi_alias(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) + sin(time), implicit_dims=ct) + eq1 = Eq(u.forward, u.forward + 1, implicit_dims=ct2) + eq2 = Eq(u.forward, u.forward + cos(time) - sin(time), implicit_dims=ct) + + op = Operator([eq0, eq1, eq2]) + op(time=5) + + cond = FindNodes(Conditional).visit(op) + assert len(cond) == 3 + # Each guard should have its own aliases for cos(time) and sin(time) + assert 'const float r0 = sin(time) + cos(time)' in str(body0(op)) + assert 'const float r1 = cos(time);' in str(body0(op)) + scalars = [i for i in FindSymbols().visit(op) if isinstance(i, Temp)] + assert len(scalars) == 3 + + def test_multi_cond_no_split(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) + # Not hoisting the init would have this equation split the time + # loop to initialize the alias for sin(time) + eq1 = Eq(u.forward, u.forward + sin(time), implicit_dims=ct2) + eq2 = Eq(u.forward, u.forward - sin(time), implicit_dims=ct) + + op = Operator([eq0, eq1, eq2]) + op(time=5) + + assert_structure( + op, + ['t', 't,x,y', 't,x,y', 't,x,y'], + 'txyxyxy' + ) + + scalars = [i for i in FindSymbols().visit(op) if isinstance(i, Temp)] + assert len(scalars) == 3 + + def test_alias_with_conditional(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(ct), implicit_dims=ct) + eq1 = Eq(u.forward, u.forward + 1, implicit_dims=ct2) + eq2 = Eq(u.forward, u.forward + cos(ct), implicit_dims=ct) + + op = Operator([eq0, eq1, eq2]) + op(time=5) + + cond = FindNodes(Conditional).visit(op) + assert len(cond) == 3 + + # Each guard should have its own alias for cos(time/ctf) + scalars = [i for i in FindSymbols().visit(op) if isinstance(i, Temp)] + assert len(scalars) == 2 + + def test_scalar_alias_interp(self): + grid = Grid(shape=(11, 11)) + time = grid.time_dim + + t_sub = ConditionalDimension(name='t_sub', parent=time, factor=3) + + f = TimeFunction(name='f', grid=grid, space_order=4) + s = SparseTimeFunction(name='s', grid=grid, npoint=1, nt=100) + + eq = Eq(f.forward, f.laplace + .002) + + rec = s.interpolate(expr=f, implicit_dims=t_sub) + + op = Operator(rec + [eq]) + + op.apply(time_M=3) + + assert np.isclose(norm(f), 254292.75, atol=0, rtol=1e-5) + assert np.isclose(norm(s), 191.44644, atol=0, rtol=1e-4) + + def test_scalar_with_cond_access(self): + grid = Grid((11, 11)) + time = grid.time_dim + + u = TimeFunction(name='u', grid=grid, time_order=2, space_order=2) + + ct = ConditionalDimension(name='ct3', parent=time, condition=Ge(time, 2)) + ct2 = ConditionalDimension(name='ct2', parent=time, factor=4) + + f1 = TimeFunction(name='f1', grid=grid, save=10, time_order=0, + dimensions=(ct,), time_dim=ct, shape=(10,)) + f1.data[:] = np.arange(10) + + eq0 = Eq(u.forward, u + cos(f1)) + eq1 = Eq(u.forward, u.forward + sin(time), implicit_dims=ct2) + eq2 = Eq(u.forward, u.forward - sin(f1)) + + op = Operator([eq0, eq1, eq2]) + + cond = FindNodes(Conditional).visit(op) + assert len(cond) == 3 + + # # Each guard should have its own alias for cos/sin(f1[time-2]) + scalars = [i for i in FindSymbols().visit(op) if isinstance(i, Temp)] + assert len(scalars) == 3 + + assert_structure( + op, + ['t', 't,x,y', 't,x,y', 't,x,y'], + 'txyxyxy' + ) + + # This would segfault without the right placement of the alias + op.apply(time_M=12) + + def test_scalar_with_cond_tinvariant(self): + grid = Grid((10, 10)) + time = grid.time_dim + dt = time.spacing + + u = TimeFunction(name='u', grid=grid, time_order=2, space_order=2) + + ct = ConditionalDimension(name='ct', parent=time, factor=2) + + eq0 = Eq(u.forward, u / dt + 1) + eq1 = Eq(u.forward, u.forward + 1/dt**2, implicit_dims=ct) + + op = Operator([eq0, eq1]) + op(time=5, dt=1) + + cond = FindNodes(Conditional).visit(op) + assert len(cond) == 1 + # One for each 1/dt 1/dt**2 + scalars = [i for i in FindSymbols().visit(op) if isinstance(i, Temp)] + assert len(scalars) == 2 + + assert_structure( + op, + ['t,x,y', 't', 't,x,y'], + 'txyxy' + ) + + # Both aliases should be hoisted outside the time loop + assert str(body0(op).body[0]) == 'const float r0 = 1.0F/dt;' + assert not body0(op).body[0].ispace + assert str(body0(op).body[1]) == 'const float r1 = 1.0F/(dt*dt);' + assert not body0(op).body[1].ispace + class TestIsoAcoustic: diff --git a/tests/test_mpi.py b/tests/test_mpi.py index 9159071997..ba3636c682 100644 --- a/tests/test_mpi.py +++ b/tests/test_mpi.py @@ -4,7 +4,7 @@ import pytest from test_dse import TestTTI -from conftest import _R, assert_blocking, assert_structure +from conftest import _R, assert_blocking, assert_structure, body0 from devito import ( NODE, Buffer, ConditionalDimension, Constant, CustomDimension, DefaultDimension, Dimension, Eq, Function, Grid, Inc, Ne, Operator, PrecomputedSparseFunction, @@ -24,12 +24,6 @@ from examples.seismic.acoustic import acoustic_setup -# Main body in Operator IET, depending on ISA -def body0(op): - bidx = 0 if 'sse' not in configuration['platform'].known_isas else 1 - return op.body.body[bidx] - - class TestDistributor: @pytest.mark.parallel(mode=[2, 4]) diff --git a/tests/test_symbolics.py b/tests/test_symbolics.py index efd55a539f..6f1675977e 100644 --- a/tests/test_symbolics.py +++ b/tests/test_symbolics.py @@ -58,8 +58,8 @@ def test_func_of_indices(): @pytest.mark.parametrize('dtype,expected', [ - (np.float32, "float r0 = 1.0F/h_x;"), - (np.float64, "double r0 = 1.0/h_x;") + (np.float32, "const float r0 = 1.0F/h_x;"), + (np.float64, "const double r0 = 1.0/h_x;") ]) def test_floatification_issue_1627(dtype, expected): """