Skip to content

Commit

Permalink
api: always make subsampling factor symbolic
Browse files Browse the repository at this point in the history
  • Loading branch information
mloubout committed Nov 6, 2023
1 parent 7c155e0 commit 7f869da
Show file tree
Hide file tree
Showing 5 changed files with 51 additions and 11 deletions.
3 changes: 2 additions & 1 deletion devito/symbolics/extended_sympy.py
Original file line number Diff line number Diff line change
Expand Up @@ -86,7 +86,8 @@ def __new__(cls, lhs, rhs, params=None):
# Perhaps it's a symbolic RHS -- but we wanna be sure it's of type int
if not hasattr(rhs, 'dtype'):
raise ValueError("Symbolic RHS `%s` lacks dtype" % rhs)
if not issubclass(rhs.dtype, np.integer):
if not issubclass(rhs.dtype, np.integer) or \
not (rhs.is_Constant and issubclass(rhs.dtype, np.integer)):
raise ValueError("Symbolic RHS `%s` must be of type `int`, found "
"`%s` instead" % (rhs, rhs.dtype))
rhs = sympify(rhs)
Expand Down
4 changes: 2 additions & 2 deletions devito/types/dense.py
Original file line number Diff line number Diff line change
Expand Up @@ -793,8 +793,8 @@ def _arg_defaults(self, alias=None):
args = ReducerMap({key.name: self._data_buffer})

# Collect default dimension arguments from all indices
for i, s in zip(self.dimensions, self.shape):
args.update(i._arg_defaults(_min=0, size=s))
for a, i, s in zip(key.dimensions, self.dimensions, self.shape):
args.update(i._arg_defaults(_min=0, size=s, alias=a))

return args

Expand Down
15 changes: 12 additions & 3 deletions devito/types/dimension.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,7 @@
from devito.tools import Pickable, is_integer
from devito.types.args import ArgProvider
from devito.types.basic import Symbol, DataSymbol, Scalar
from devito.types.constant import Constant

__all__ = ['Dimension', 'SpaceDimension', 'TimeDimension', 'DefaultDimension',
'CustomDimension', 'SteppingDimension', 'SubDimension', 'ConditionalDimension',
Expand Down Expand Up @@ -848,13 +849,21 @@ def __init_finalize__(self, name, parent=None, factor=None, condition=None,

super().__init_finalize__(name, parent)

self._factor = factor
# Always make the factor symbolic to allow overrides with different factor.
if factor is None:
self._factor = None
elif is_integer(factor):
self._factor = Constant(name="%sf" % name, value=factor, dtype=np.int32)
elif factor.is_Constant and is_integer(factor.data):
self._factor = factor
else:
raise ValueError("factor must be an integer or integer Constant")
self._condition = condition
self._indirect = indirect

@property
def spacing(self):
return self.factor * self.parent.spacing
return self.factor.data * self.parent.spacing

@property
def factor(self):
Expand Down Expand Up @@ -890,7 +899,7 @@ def _arg_defaults(self, _min=None, size=None, alias=None):
return defaults
try:
# Is it a symbolic factor?
factor = defaults[dim._factor.name] = dim._factor.data
factor = defaults[dim._factor.name] = self._factor.data
except AttributeError:
factor = dim._factor

Expand Down
2 changes: 1 addition & 1 deletion devito/types/grid.py
Original file line number Diff line number Diff line change
Expand Up @@ -280,7 +280,7 @@ def spacing_map(self):
# Special case subsampling: `Grid.dimensions` -> (xb, yb, zb)`
# where `xb, yb, zb` are ConditionalDimensions whose parents
# are SpaceDimensions
mapper[d.root.spacing] = s/self.dtype(d.factor)
mapper[d.root.spacing] = s/self.dtype(d.factor.data)
elif d.is_Space:
# Typical case: `Grid.dimensions` -> (x, y, z)` where `x, y, z` are
# the SpaceDimensions
Expand Down
38 changes: 34 additions & 4 deletions tests/test_dimension.py
Original file line number Diff line number Diff line change
Expand Up @@ -793,6 +793,34 @@ def test_overrides(self):
assert np.all([np.allclose(uk.data[i], i*fk)
for i in range((nt+fk-1)//fk)])

def test_overrides_newfact(self):
nt = 19
grid = Grid(shape=(11, 11))
time = grid.time_dim

u = TimeFunction(name='u', grid=grid)
assert(grid.stepping_dim in u.indices)

f1, f2 = 4, 5
n1, n2 = (nt+f1-1)//f1, (nt+f2-1)//f2
t1 = ConditionalDimension('t_sub1', parent=time, factor=f1)
t2 = ConditionalDimension('t_sub2', parent=time, factor=f2)
u1 = TimeFunction(name='usave1', grid=grid, save=n1, time_dim=t1)
u2 = TimeFunction(name='usave2', grid=grid, save=n2, time_dim=t2)
assert(t1 in u1.indices)
assert(t2 in u2.indices)

eqns = [Eq(u.forward, u + 1.), Eq(u1, u)]
op = Operator(eqns)
op.apply(u=u, usave1=u1, time_M=nt-2)
u.data.fill(0)
op.apply(u=u, usave1=u2, time_M=nt-2)

assert np.all(np.allclose(u.data[(nt-1) % 3], nt-1))
for (uk, fk) in zip((u1, u2), (f1, f2)):
assert np.all([np.allclose(uk.data[i], i*fk)
for i in range((nt+fk-1)//fk)])

def test_basic_shuffles(self):
"""
Like ``test_basic``, but with different equation orderings. Nevertheless,
Expand Down Expand Up @@ -871,10 +899,10 @@ def test_time_subsampling_fd(self):
usave = TimeFunction(name='usave', grid=grid, save=(nt+factor-1)//factor,
time_dim=time_subsampled, time_order=2)

dx2 = [indexify(i) for i in retrieve_functions(usave.dt2.evaluate)]
assert dx2 == [usave[time_subsampled - 1, x, y],
dx2 = {indexify(i) for i in retrieve_functions(usave.dt2.evaluate)}
assert dx2 == {usave[time_subsampled - 1, x, y],
usave[time_subsampled + 1, x, y],
usave[time_subsampled, x, y]]
usave[time_subsampled, x, y]}

def test_issue_1592(self):
grid = Grid(shape=(11, 11))
Expand All @@ -886,7 +914,9 @@ def test_issue_1592(self):
op = Operator(Eq(v.forward, v.dx))
op.apply(time=6)
exprs = FindNodes(Expression).visit(op)
assert exprs[-1].expr.lhs.indices[0] == IntDiv(time, 2) + 1
assert exprs[-1].expr.lhs.indices[0] == IntDiv(time, time_sub.factor) + 1
assert time_sub.factor.data == 2
assert time_sub.factor.is_Constant

def test_issue_1753(self):
grid = Grid(shape=(3, 3, 3))
Expand Down

0 comments on commit 7f869da

Please sign in to comment.