Merge branch 'master' into pull-request-cleanup

devito/__init__.py

@@ -1,14 +1,12 @@
 from __future__ import absolute_import
 
 from devito.base import *  # noqa
-from devito.finite_difference import *  # noqa
 from devito.dimension import *  # noqa
+from devito.equation import *  # noqa
+from devito.finite_difference import *  # noqa
 from devito.grid import *  # noqa
-from devito.function import Forward, Backward  # noqa
 from devito.logger import error, warning, info, set_log_level  # noqa
 from devito.parameters import *  # noqa
-from devito.region import *  # noqa
-from devito.symbolics import *  # noqa
 from devito.tools import *  # noqa
 
 from devito.compiler import compiler_registry

devito/dimension.py

@@ -27,13 +27,11 @@ class Dimension(AbstractSymbol):
     potential iteration space.
 
     :param name: Name of the dimension symbol.
-    :param reverse: Optional, Traverse dimension in reverse order (default False)
     :param spacing: Optional, symbol for the spacing along this dimension.
     """
 
     def __new__(cls, name, **kwargs):
         newobj = sympy.Symbol.__new__(cls, name)
-        newobj._reverse = kwargs.get('reverse', False)
         newobj._spacing = kwargs.get('spacing', Scalar(name='h_%s' % name))
         return newobj
 
@@ -86,27 +84,16 @@ def start_name(self):
     def end_name(self):
         return "%s_e" % self.name
 
-    @property
-    def reverse(self):
-        return self._reverse
-
     @property
     def spacing(self):
         return self._spacing
 
-    @reverse.setter
-    def reverse(self, val):
-        # TODO: this is an outrageous hack. TimeFunctions are updating this value
-        # at construction time.
-        self._reverse = val
-
     @property
     def base(self):
         return self
 
     def _hashable_content(self):
-        return super(Dimension, self)._hashable_content() +\
-            (self.reverse, self.spacing)
+        return super(Dimension, self)._hashable_content() + (self.spacing,)
 
     def argument_defaults(self, size=None):
         """
@@ -148,7 +135,6 @@ class SpaceDimension(Dimension):
     symbols.
 
     :param name: Name of the dimension symbol.
-    :param reverse: Traverse dimension in reverse order (default False)
     :param spacing: Optional, symbol for the spacing along this dimension.
     """
 
@@ -163,7 +149,6 @@ class TimeDimension(Dimension):
     time dimensions should inherit from :class:`TimeDimension`.
 
     :param name: Name of the dimension symbol.
-    :param reverse: Traverse dimension in reverse order (default False)
     :param spacing: Optional, symbol for the spacing along this dimension.
     """
 
@@ -193,10 +178,6 @@ def __new__(cls, name, parent, **kwargs):
     def parent(self):
         return self._parent
 
-    @property
-    def reverse(self):
-        return self.parent.reverse
-
     @property
     def spacing(self):
         return self.parent.spacing

devito/dle/blocking_utils.py

@@ -201,11 +201,9 @@ class IterationFold(Iteration):
 
     is_IterationFold = True
 
-    def __init__(self, nodes, dimension, limits, index=None, offsets=None,
-                 properties=None, pragmas=None, uindices=None, folds=None):
-        super(IterationFold, self).__init__(nodes, dimension, limits, index, offsets,
-                                            properties, uindices, pragmas)
-        self.folds = folds
+    def __init__(self, *args, **kwargs):
+        self.folds = kwargs.pop('folds', None)
+        super(IterationFold, self).__init__(*args, **kwargs)
 
     def __repr__(self):
         properties = ""

devito/dse/backends/advanced.py

@@ -2,7 +2,7 @@
 
 from collections import OrderedDict
 
-from devito.ir import Cluster, ClusterGroup, groupby
+from devito.ir import Cluster, ClusterGroup, IterationSpace, groupby
 from devito.dse.aliases import collect
 from devito.dse.backends import BasicRewriter, dse_pass
 from devito.symbolics import Eq, estimate_cost, xreplace_constrained, iq_timeinvariant
@@ -140,9 +140,12 @@ def _eliminate_inter_stencil_redundancies(self, cluster, template, **kwargs):
             function = make(c)
             # Build new Cluster
             expression = Eq(Indexed(function, *indices), origin)
-            ispace = cluster.ispace.subtract(alias.anti_stencil.boxify().negate())
+            intervals, sub_iterators, directions = cluster.ispace.args
+            # Adjust intervals
+            intervals = intervals.subtract(alias.anti_stencil.boxify().negate())
             if all(time_invariants[i] for i in alias.aliased):
-                ispace = ispace.drop([i.dim for i in ispace.intervals if i.dim.is_Time])
+                intervals = intervals.drop([i for i in intervals.dimensions if i.is_Time])
+            ispace = IterationSpace(intervals, sub_iterators, directions)
             alias_clusters.append(Cluster([expression], ispace))
             # Update substitution rules
             for aliased, distance in alias.with_distance:

devito/equation.py

@@ -0,0 +1,62 @@
+"""User API to specify equations."""
+
+import sympy
+
+__all__ = ['Eq', 'Inc', 'DOMAIN', 'INTERIOR']
+
+
+class Eq(sympy.Eq):
+
+    """
+    A :class:`sympy.Eq` that accepts the additional keyword parameter ``region``.
+
+    The ``region``, an object of type :class:`Region`, may be used to restrict
+    the execution of the equation to a sub-domain.
+    """
+
+    is_Increment = False
+
+    def __new__(cls, *args, **kwargs):
+        kwargs['evaluate'] = False
+        region = kwargs.pop('region', DOMAIN)
+        obj = sympy.Eq.__new__(cls, *args, **kwargs)
+        obj._region = region
+        return obj
+
+
+class Inc(Eq):
+
+    """
+    A :class:`Eq` performing a linear increment.
+    """
+
+    is_Increment = True
+
+
+class Region(object):
+
+    """
+    A region of the computational domain over which a :class:`Function` is
+    discretized.
+    """
+
+    def __init__(self, name):
+        self._name = name
+
+    def __repr__(self):
+        return self._name
+
+    def __eq__(self, other):
+        return isinstance(other, Region) and self._name == other._name
+
+
+DOMAIN = Region('DOMAIN')
+"""
+Represent the physical domain of the PDE; that is, domain = boundary + interior
+"""
+
+INTERIOR = Region('INTERIOR')
+"""
+Represent the physical interior domain of the PDE; that is, PDE boundaries are
+not included.
+"""

devito/function.py

@@ -6,44 +6,23 @@
 import numpy as np
 from psutil import virtual_memory
 
-from devito.parameters import configuration
-from devito.logger import debug, error, warning
-from devito.data import Data, first_touch
+from devito.arguments import ArgumentMap
 from devito.cgen_utils import INT, FLOAT
+from devito.data import Data, first_touch
 from devito.dimension import Dimension
-from devito.types import SymbolicFunction, AbstractCachedSymbol
-from devito.arguments import ArgumentMap
+from devito.equation import Eq, Inc
 from devito.finite_difference import (centered, cross_derivative,
                                       first_derivative, left, right,
                                       second_derivative, generic_derivative,
                                       second_cross_derivative)
-from devito.symbolics import Eq, Inc, indexify, retrieve_indexed
+from devito.logger import debug, error, warning
+from devito.parameters import configuration
+from devito.symbolics import indexify, retrieve_indexed
+from devito.types import SymbolicFunction, AbstractCachedSymbol
 from devito.tools import EnrichedTuple
 
 __all__ = ['Constant', 'Function', 'TimeFunction', 'SparseFunction',
-           'SparseTimeFunction', 'Forward', 'Backward']
-
-
-class TimeAxis(object):
-    """Direction in which to advance the time index on
-    :class:`TimeFunction` objects.
-
-    :param axis: Either 'Forward' or 'Backward'
-    """
-
-    def __init__(self, axis):
-        assert axis in ['Forward', 'Backward']
-        self._axis = {'Forward': 1, 'Backward': -1}[axis]
-
-    def __eq__(self, other):
-        return self._axis == other._axis
-
-    def __repr__(self):
-        return {-1: 'Backward', 1: 'Forward'}[self._axis]
-
-
-Forward = TimeAxis('Forward')
-Backward = TimeAxis('Backward')
+           'SparseTimeFunction']
 
 
 class Constant(AbstractCachedSymbol):
@@ -333,6 +312,7 @@ def argument_defaults(self, alias=None):
         # Collect default dimension arguments from all indices
         for i, s, o in zip(self.indices, self.shape, self.staggered):
             args.update(i.argument_defaults(size=s+o))
+
         return args
 
     def argument_values(self, alias=None, **kwargs):
@@ -348,18 +328,20 @@ def argument_values(self, alias=None, **kwargs):
         # Add value override for own data if it is provided
         if self.name in kwargs:
             new = kwargs.pop(self.name)
+            if len(new.shape) != self.ndim:
+                raise ValueError("Array shape %s does not match" % (new.shape, ) +
+                                 "dimensions %s" % (self.indices, ))
             if isinstance(new, TensorFunction):
                 # Set new values and re-derive defaults
                 values[key] = new._data_buffer
                 values.update(new.argument_defaults(alias=key).reduce_all())
             else:
                 # We've been provided a pure-data replacement (array)
                 values[key] = new
-                # TODO: Re-derive defaults from shape of array
+                # Add value overrides for all associated dimensions
+                for i, s, o in zip(self.indices, new.shape, self.staggered):
+                    values.update(i.argument_defaults(size=s+o))
 
-        # Add value overrides for all associated dimensions
-        for i in self.indices:
-            values.update(i.argument_values(**kwargs))
         return values
 
 
@@ -659,6 +641,7 @@ def __init__(self, *args, **kwargs):
 
             self.time_dim = kwargs.get('time_dim')
             self.time_order = kwargs.get('time_order', 1)
+            self.save = type(kwargs.get('save', None) or None)
             if not isinstance(self.time_order, int):
                 raise ValueError("'time_order' must be int")
 
@@ -744,6 +727,23 @@ def dt2(self):
 
         return self.diff(_t, _t).as_finite_difference(indt)
 
+    def argument_values(self, alias=None, **kwargs):
+        """
+        Returns a map of argument values after evaluating user input.
+
+        :param kwargs: Dictionary of user-provided argument overrides.
+        :param alias: (Optional) name under which to store values.
+        """
+        # Check if data has the right dimension
+        if self.name in kwargs:
+            new = kwargs.get(self.name)
+            if isinstance(new, TimeFunction) and new.save != self.save:
+                raise TypeError("Incorrect value encountered, save should be %s" %
+                                self.save)
+
+        values = super(TimeFunction, self).argument_values(alias=alias, **kwargs)
+        return values
+
 
 class SparseFunction(TensorFunction):
     """
@@ -816,7 +816,7 @@ def __indices_setup__(cls, **kwargs):
         if dimensions is not None:
             return dimensions
         else:
-            return (Dimension(name='p'),)
+            return (Dimension(name='p_%s' % kwargs["name"]),)
 
     @classmethod
     def __shape_setup__(cls, **kwargs):
@@ -1022,16 +1022,20 @@ def argument_values(self, alias=None, **kwargs):
         :param alias: (Optional) name under which to store values.
         """
         # Take a copy of the replacement before super pops it from kwargs
+
         new = kwargs.get(self.name)
-        values = super(SparseFunction, self).argument_values(alias=alias, **kwargs)
+        key = alias or self.name
+
+        values = super(SparseFunction, self).argument_values(alias=key, **kwargs)
+
         if new is not None and isinstance(new, SparseFunction):
             # If we've been replaced with a SparseFunction,
             # we need to re-derive defaults and values...
-            values.update(new.argument_defaults(alias=alias).reduce_all())
+            values.update(new.argument_defaults(alias=key).reduce_all())
             values.update(new.coordinates.argument_defaults(alias=self.coordinates.name))
         else:
             # ..., but if not, we simply need to recurse over children.
-            values.update(self.coordinates.argument_values(alias=alias, **kwargs))
+            values.update(self.coordinates.argument_values(alias=key, **kwargs))
         return values
 
 
@@ -1080,7 +1084,7 @@ def __indices_setup__(cls, **kwargs):
         if dimensions is not None:
             return dimensions
         else:
-            return (kwargs['grid'].time_dim, Dimension(name='p'))
+            return (kwargs['grid'].time_dim, Dimension(name='p_%s' % kwargs["name"]))
 
     @classmethod
     def __shape_setup__(cls, **kwargs):