[Yt-svn] commit/yt: MatthewTurk: A few more changes for speed in smoothed covering grids. Fixing the previous

Bitbucket commits-noreply at bitbucket.org
Thu Sep 8 11:07:36 PDT 2011


1 new changeset in yt:

http://bitbucket.org/yt_analysis/yt/changeset/fe04156a8573/
changeset:   fe04156a8573
branch:      yt
user:        MatthewTurk
date:        2011-09-08 20:07:26
summary:     A few more changes for speed in smoothed covering grids.  Fixing the previous
regression brought it back to a net 2x speedup (but, accurate again) and this
shaves about 20% off that.
affected #:  2 files (513 bytes)

--- a/yt/data_objects/data_containers.py	Thu Sep 08 11:13:51 2011 -0600
+++ b/yt/data_objects/data_containers.py	Thu Sep 08 14:07:26 2011 -0400
@@ -3104,6 +3104,9 @@
     _type_name = "smoothed_covering_grid"
     @wraps(AMRCoveringGridBase.__init__)
     def __init__(self, *args, **kwargs):
+        self._base_dx = (
+              (self.pf.domain_right_edge - self.pf.domain_left_edge) /
+               self.pf.domain_dimensions.astype("float64"))
         AMRCoveringGridBase.__init__(self, *args, **kwargs)
         self._final_start_index = self.global_startindex
 
@@ -3155,10 +3158,10 @@
             if self._use_pbar: pbar.finish()
 
     def _update_level_state(self, level, field = None):
-        dx = ((self.pf.domain_right_edge - self.pf.domain_left_edge) /
-               self.pf.domain_dimensions.astype("float64"))
-        dx /= self.pf.refine_by**level
-        for ax, v in zip('xyz', dx): self['cd%s'%ax] = v
+        dx = self._base_dx / self.pf.refine_by**level
+        self.data['cdx'] = dx[0]
+        self.data['cdy'] = dx[1]
+        self.data['cdz'] = dx[2]
         LL = self.left_edge - self.pf.domain_left_edge
         self._old_global_startindex = self.global_startindex
         self.global_startindex = na.rint(LL / dx).astype('int64') - 1
@@ -3167,37 +3170,28 @@
         if level == 0 and self.level > 0:
             # We use one grid cell at LEAST, plus one buffer on all sides
             idims = na.rint((self.right_edge-self.left_edge)/dx).astype('int64') + 2
-            self[field] = na.zeros(idims,dtype='float64')-999
+            self.data[field] = na.zeros(idims,dtype='float64')-999
             self._cur_dims = idims.astype("int32")
         elif level == 0 and self.level == 0:
             DLE = self.pf.domain_left_edge
             self.global_startindex = na.array(na.floor(LL/ dx), dtype='int64')
             idims = na.rint((self.right_edge-self.left_edge)/dx).astype('int64')
-            self[field] = na.zeros(idims,dtype='float64')-999
+            self.data[field] = na.zeros(idims,dtype='float64')-999
             self._cur_dims = idims.astype("int32")
 
     def _refine(self, dlevel, field):
         rf = float(self.pf.refine_by**dlevel)
 
-        old_dims = na.array(self[field].shape) - 1
-        old_left = (self._old_global_startindex + 0.5) * rf 
-        old_right = rf*old_dims + old_left
-        dx = na.array([self['cd%s' % ax] for ax in 'xyz'], dtype='float64')
-        new_dims = na.rint((self.right_edge-self.left_edge)/dx).astype('int64') + 2
+        input_left = (self._old_global_startindex + 0.5) * rf 
+        dx = na.fromiter((self['cd%s' % ax] for ax in 'xyz'), count=3, dtype='float64')
+        output_dims = na.rint((self.right_edge-self.left_edge)/dx).astype('int32') + 2
 
-        self._cur_dims = new_dims.astype("int32")
+        self._cur_dims = output_dims
 
-        output_field = na.zeros(new_dims, dtype="float64")
-        input_bounds = na.array([[old_left[0], old_right[0]],
-                                 [old_left[1], old_right[1]],
-                                 [old_left[2], old_right[2]]])
-        new_left = self.global_startindex + 0.5
-        new_right = new_left + new_dims - 1
-        output_bounds = na.array([[new_left[0], new_right[0]],
-                                  [new_left[1], new_right[1]],
-                                  [new_left[2], new_right[2]]])
-        ghost_zone_interpolate(self[field], input_bounds,
-                               output_field, output_bounds)
+        output_field = na.zeros(output_dims, dtype="float64")
+        output_left = self.global_startindex + 0.5
+        ghost_zone_interpolate(rf, self[field], input_left,
+                               output_field, output_left)
         self[field] = output_field
 
     def _get_data_from_grid(self, grid, fields, level):


--- a/yt/utilities/_amr_utils/Interpolators.pyx	Thu Sep 08 11:13:51 2011 -0600
+++ b/yt/utilities/_amr_utils/Interpolators.pyx	Thu Sep 08 14:07:26 2011 -0400
@@ -124,35 +124,38 @@
 @cython.cdivision(True)
 @cython.wraparound(False)
 @cython.boundscheck(False)
-def ghost_zone_interpolate(np.ndarray[np.float64_t, ndim=3] input_field,
-                           np.ndarray[np.float64_t, ndim=2] input_bounds,
+def ghost_zone_interpolate(int rf,
+                           np.ndarray[np.float64_t, ndim=3] input_field,
+                           np.ndarray[np.float64_t, ndim=1] input_left,
                            np.ndarray[np.float64_t, ndim=3] output_field,
-                           np.ndarray[np.float64_t, ndim=2] output_bounds):
+                           np.ndarray[np.float64_t, ndim=1] output_left):
     cdef int oi, oj, ok
     cdef int ii, ij, ik
     cdef np.float64_t xp, xm, yp, ym, zp, zm
     cdef np.float64_t ods[3], ids[3], iids[3]
-    cdef np.float64_t opos[3], ropos[3]
+    cdef np.float64_t opos[3], ropos[3], temp
     cdef int i, j
     for i in range(3):
-        ids[i] = (input_bounds[i,1] - input_bounds[i,0])/(input_field.shape[i]-1)
-        ods[i] = (output_bounds[i,1] - output_bounds[i,0])/(output_field.shape[i]-1)
+        temp = input_left[i] + (rf * (input_field.shape[i] - 1))
+        ids[i] = (temp - input_left[i])/(input_field.shape[i]-1)
+        temp = output_left[i] + output_field.shape[i] - 1
+        ods[i] = (temp - output_left[i])/(output_field.shape[i]-1)
         iids[i] = 1.0/ids[i]
-    opos[0] = output_bounds[0,0]
+    opos[0] = output_left[0]
     for oi in range(output_field.shape[0]):
-        ropos[0] = ((opos[0] - input_bounds[0,0]) * iids[0])
+        ropos[0] = ((opos[0] - input_left[0]) * iids[0])
         ii = iclip(<int> ropos[0], 0, input_field.shape[0] - 2)
         xp = ropos[0] - ii
         xm = 1.0 - xp
-        opos[1] = output_bounds[1,0]
+        opos[1] = output_left[1]
         for oj in range(output_field.shape[1]):
-            ropos[1] = ((opos[1] - input_bounds[1,0]) * iids[1])
+            ropos[1] = ((opos[1] - input_left[1]) * iids[1])
             ij = iclip(<int> ropos[1], 0, input_field.shape[1] - 2)
             yp = ropos[1] - ij
             ym = 1.0 - yp
-            opos[2] = output_bounds[2,0]
+            opos[2] = output_left[2]
             for ok in range(output_field.shape[2]):
-                ropos[2] = ((opos[2] - input_bounds[2,0]) * iids[2])
+                ropos[2] = ((opos[2] - input_left[2]) * iids[2])
                 ik = iclip(<int> ropos[2], 0, input_field.shape[2] - 2)
                 zp = ropos[2] - ik
                 zm = 1.0 - zp

Repository URL: https://bitbucket.org/yt_analysis/yt/

--

This is a commit notification from bitbucket.org. You are receiving
this because you have the service enabled, addressing the recipient of
this email.



More information about the yt-svn mailing list