[Yt-svn] commit/yt: 2 new changesets
Bitbucket
commits-noreply at bitbucket.org
Tue Jun 7 17:23:33 PDT 2011
2 new changesets in yt:
http://bitbucket.org/yt_analysis/yt/changeset/c45d638fb307/
changeset: c45d638fb307
branch: yt
user: MatthewTurk
date: 2011-06-08 01:03:15
summary: Keep all trees resident in memory in Ramses to save time at the expense of
memory. (For now.)
affected #: 2 files (268 bytes)
--- a/yt/frontends/ramses/_ramses_reader.pyx Tue Jun 07 16:00:46 2011 -0700
+++ b/yt/frontends/ramses/_ramses_reader.pyx Tue Jun 07 19:03:15 2011 -0400
@@ -397,10 +397,10 @@
# We now have to get our field names to fill our array
self.trees = <RAMSES_tree**>\
malloc(sizeof(RAMSES_tree*) * self.rsnap.m_header.ncpu)
+ for ii in range(self.ndomains): self.trees[ii] = NULL
self.hydro_datas = <RAMSES_hydro_data ***>\
malloc(sizeof(RAMSES_hydro_data**) * self.rsnap.m_header.ncpu)
self.ndomains = self.rsnap.m_header.ncpu
- #for ii in range(self.ndomains): self.trees[ii] = NULL
# Note we don't do ncpu + 1
for idomain in range(self.rsnap.m_header.ncpu):
# we don't delete local_tree
@@ -415,8 +415,8 @@
new RAMSES_hydro_data(deref(local_tree))
self.trees[idomain] = local_tree
# We do not delete the final snapshot, which we'll use later
- if idomain + 1 < self.rsnap.m_header.ncpu:
- del local_hydro_data
+ #if idomain + 1 < self.rsnap.m_header.ncpu:
+ # del local_hydro_data
# Only once, we read all the field names
self.nfields = local_hydro_data.m_nvars
cdef string *field_name
@@ -434,7 +434,6 @@
self.field_names.append(field_name.c_str())
self.field_ind[self.field_names[-1]] = ifield
# This all needs to be cleaned up in the deallocator
- del local_hydro_data
def __dealloc__(self):
import traceback; traceback.print_stack()
@@ -473,21 +472,20 @@
cdef np.ndarray[np.int64_t, ndim=1] cell_count
cell_count = np.zeros(self.rsnap.m_header.levelmax + 1, 'int64')
cdef int local_count = 0
+ cdef int tree_count
for idomain in range(1, self.rsnap.m_header.ncpu + 1):
- local_tree = new RAMSES_tree(deref(self.rsnap), idomain,
- self.rsnap.m_header.levelmax, 0)
- local_tree.read()
- local_hydro_data = new RAMSES_hydro_data(deref(local_tree))
+ local_tree = self.trees[idomain - 1]
for ilevel in range(local_tree.m_maxlevel + 1):
local_count = 0
+ tree_count = 0
local_level = &local_tree.m_AMR_levels[ilevel]
grid_it = local_tree.begin(ilevel)
grid_end = local_tree.end(ilevel)
while grid_it != grid_end:
local_count += (grid_it.get_domain() == idomain)
+ tree_count += 1
grid_it.next()
cell_count[ilevel] += local_count
- del local_tree, local_hydro_data
return cell_count
@@ -582,10 +580,7 @@
cdef np.ndarray[np.int64_t, ndim=1] level_cell_counts
level_cell_counts = np.zeros(self.rsnap.m_header.levelmax + 1, 'int64')
for idomain in range(1, self.rsnap.m_header.ncpu + 1):
- local_tree = new RAMSES_tree(deref(self.rsnap), idomain,
- self.rsnap.m_header.levelmax, 0)
- local_tree.read()
- local_hydro_data = new RAMSES_hydro_data(deref(local_tree))
+ local_tree = self.trees[idomain - 1]
for ilevel in range(local_tree.m_maxlevel + 1):
# this gets overwritten for every domain, which is okay
level_cell_counts[ilevel] = grid_ind
@@ -626,7 +621,6 @@
grid_ind += 1
grid_aind += 1
grid_it.next()
- del local_tree, local_hydro_data
def read_oct_grid(self, char *field, int level, int domain, int grid_id):
@@ -715,13 +709,29 @@
to_fill += 1
return to_fill
+#def recursive_patch_splitting(ProtoSubgrid psg,
+# np.ndarray[np.int64_t, ndim=1] dims,
+# np.ndarray[np.int64_t, ndim=1] inds,
+# np.ndarray[np.int64_t, ndim=2] left_index,
+# np.ndarray[np.int64_t, ndim=2] right_index,
+# np.ndarray[np.int64_t, ndim=2] gdims,
+# np.ndarray[np.int64_t, ndim=2] fl,
+# int num_deep = 0):
+# cdef float min_eff = 0.1
+# if num_deep > 40:
+# psg.efficiency = min_eff
+# return [psg]
+# if psg.efficiency > min_eff or psg.efficiency < 0.0:
+# return [psg]
+# cdef
+#
cdef class ProtoSubgrid:
cdef np.int64_t *signature[3]
cdef np.int64_t left_edge[3]
cdef np.int64_t right_edge[3]
cdef np.int64_t dimensions[3]
cdef public np.float64_t efficiency
- cdef public object sigs
+ cdef np.int64_t *sigs[3]
cdef public object grid_file_locations
cdef public object dd
@@ -740,7 +750,6 @@
cdef int l0, r0, l1, r1, l2, r2, i0, i1, i2
cdef np.int64_t temp_l[3], temp_r[3], ncells
cdef np.float64_t efficiency
- self.sigs = []
for i in range(3):
temp_l[i] = left_index[i] + dimensions[i]
temp_r[i] = left_index[i]
@@ -764,13 +773,14 @@
if self.dimensions[i] <= 0:
self.efficiency = -1.0
return
- self.sigs.append(np.zeros(self.dimensions[i], 'int64'))
- #print self.sigs[0].size, self.sigs[1].size, self.sigs[2].size
+ self.sigs[i] = <np.int64_t *> malloc(
+ sizeof(np.int64_t) * self.dimensions[i])
+ for gi in range(self.dimensions[i]): self.sigs[i][gi] = 0
# My guess is that this whole loop could be done more efficiently.
# However, this is clear and straightforward, so it is a good first
# pass.
- cdef np.ndarray[np.int64_t, ndim=1] sig0, sig1, sig2
+ cdef np.int64_t *sig0, *sig1, *sig2
sig0 = self.sigs[0]
sig1 = self.sigs[1]
sig2 = self.sigs[2]
@@ -813,6 +823,11 @@
#print "Efficiency is %0.3e" % (efficiency)
self.efficiency = efficiency
+ def __dealloc__(self):
+ free(self.sigs[0])
+ free(self.sigs[1])
+ free(self.sigs[2])
+
@cython.boundscheck(False)
@cython.wraparound(False)
def find_split(self):
@@ -821,7 +836,7 @@
cdef np.ndarray[ndim=1, dtype=np.int64_t] axes
cdef np.int64_t strength, zcstrength, zcp
axes = np.argsort(self.dd)[::-1]
- cdef np.ndarray[np.int64_t] sig
+ cdef np.int64_t *sig
for axi in range(3):
ax = axes[axi]
center = self.dimensions[ax] / 2
--- a/yt/frontends/ramses/data_structures.py Tue Jun 07 16:00:46 2011 -0700
+++ b/yt/frontends/ramses/data_structures.py Tue Jun 07 19:03:15 2011 -0400
@@ -203,9 +203,9 @@
sums = na.zeros(3, dtype='int64')
mylog.info("Final grid count: %s", len(self.proto_grids[level]))
if len(self.proto_grids[level]) == 1: continue
- for g in self.proto_grids[level]:
- sums += [s.sum() for s in g.sigs]
- assert(na.all(sums == dims.prod(axis=1).sum()))
+ #for g in self.proto_grids[level]:
+ # sums += [s.sum() for s in g.sigs]
+ #assert(na.all(sums == dims.prod(axis=1).sum()))
self.num_grids = sum(len(l) for l in self.proto_grids)
num_deep = 0
http://bitbucket.org/yt_analysis/yt/changeset/0a029dc99f32/
changeset: 0a029dc99f32
branch: yt
user: MatthewTurk
date: 2011-06-08 02:22:57
summary: More improvements to RAMSES reader; another factor of a few speedup by removing
some array copies and pushing to Cython some of the bigger tasks.
affected #: 3 files (1.1 KB)
--- a/yt/frontends/ramses/_ramses_reader.pyx Tue Jun 07 19:03:15 2011 -0400
+++ b/yt/frontends/ramses/_ramses_reader.pyx Tue Jun 07 20:22:57 2011 -0400
@@ -741,8 +741,6 @@
np.ndarray[np.int64_t, ndim=1] left_index,
np.ndarray[np.int64_t, ndim=1] dimensions,
np.ndarray[np.int64_t, ndim=2] left_edges,
- np.ndarray[np.int64_t, ndim=2] right_edges,
- np.ndarray[np.int64_t, ndim=2] grid_dimensions,
np.ndarray[np.int64_t, ndim=2] grid_file_locations):
# This also includes the shrinking step.
cdef int i, ci, ng = left_edges.shape[0]
@@ -756,16 +754,16 @@
self.signature[i] = NULL
for gi in range(ng):
if left_edges[gi,0] > left_index[0]+dimensions[0] or \
- right_edges[gi,0] < left_index[0] or \
+ left_edges[gi,0] + 2 < left_index[0] or \
left_edges[gi,1] > left_index[1]+dimensions[1] or \
- right_edges[gi,1] < left_index[1] or \
+ left_edges[gi,1] + 2 < left_index[1] or \
left_edges[gi,2] > left_index[2]+dimensions[2] or \
- right_edges[gi,2] < left_index[2]:
+ left_edges[gi,2] + 2 < left_index[2]:
#print "Skipping grid", gi, "which lies outside out box"
continue
for i in range(3):
temp_l[i] = i64min(left_edges[gi,i], temp_l[i])
- temp_r[i] = i64max(right_edges[gi,i], temp_r[i])
+ temp_r[i] = i64max(left_edges[gi,i] + 2, temp_r[i])
for i in range(3):
self.left_edge[i] = i64max(temp_l[i], left_index[i])
self.right_edge[i] = i64min(temp_r[i], left_index[i] + dimensions[i])
@@ -790,15 +788,15 @@
for gi in range(ng):
used = 0
nnn = 0
- for l0 in range(grid_dimensions[gi, 0]):
+ for l0 in range(2):
i0 = left_edges[gi, 0] + l0
if i0 < self.left_edge[0]: continue
if i0 >= self.right_edge[0]: break
- for l1 in range(grid_dimensions[gi, 1]):
+ for l1 in range(2):
i1 = left_edges[gi, 1] + l1
if i1 < self.left_edge[1]: continue
if i1 >= self.right_edge[1]: break
- for l2 in range(grid_dimensions[gi, 2]):
+ for l2 in range(2):
i2 = left_edges[gi, 2] + l2
if i2 < self.left_edge[2]: continue
if i2 >= self.right_edge[2]: break
@@ -969,3 +967,36 @@
hilbert_indices[o] = h
return hilbert_indices
+ at cython.boundscheck(False)
+ at cython.wraparound(False)
+def get_array_indices_lists(np.ndarray[np.int64_t, ndim=1] ind,
+ np.ndarray[np.int64_t, ndim=1] uind):
+ cdef np.ndarray[np.int64_t, ndim=1] count = np.zeros(uind.shape[0], 'int64')
+ cdef int n, i
+ cdef np.int64_t mi, mui
+ for i in range(ind.shape[0]):
+ mi = ind[i]
+ for n in range(uind.shape[0]):
+ if uind[n] == mi:
+ count[n] += 1
+ break
+ cdef np.int64_t **inds
+ inds = <np.int64_t **> malloc(sizeof(np.int64_t *) * uind.shape[0])
+ cdef int *li = <int *> malloc(sizeof(int) * uind.shape[0])
+ cdef np.ndarray[np.int64_t, ndim=1] indices
+ all_indices = []
+ for n in range(uind.shape[0]):
+ indices = np.zeros(count[n], 'int64')
+ all_indices.append(indices)
+ inds[n] = <np.int64_t *> indices.data
+ li[n] = 0
+ for i in range(ind.shape[0]):
+ mi = ind[i]
+ for n in range(uind.shape[0]):
+ if uind[n] == mi:
+ inds[n][li[n]] = i
+ li[n] += 1
+ break
+ free(inds) # not inds[...]
+ free(li)
+ return all_indices
--- a/yt/frontends/ramses/data_structures.py Tue Jun 07 19:03:15 2011 -0400
+++ b/yt/frontends/ramses/data_structures.py Tue Jun 07 20:22:57 2011 -0400
@@ -38,6 +38,8 @@
from .fields import RAMSESFieldContainer
from yt.utilities.definitions import \
mpc_conversion
+from yt.utilities.amr_utils import \
+ get_box_grids_level
from yt.utilities.io_handler import \
io_registry
@@ -166,7 +168,6 @@
# left_index is integers of the index, with respect to this level
left_index = na.rint((ogrid_left_edge[ggi,:]) * nd / DW ).astype('int64')
# we've got octs, so it's +2
- right_index = left_index + 2
pbar = get_pbar("Re-gridding ", left_index.shape[0])
dlp = [None, None, None]
i = 0
@@ -180,24 +181,27 @@
# Strictly speaking, we don't care about the index of any
# individual oct at this point. So we can then split them up.
unique_indices = na.unique(hilbert_indices)
- for curve_index in unique_indices:
+ print "Level % 2i has % 10i unique indices for %0.3e octs" % (
+ level, unique_indices.size, hilbert_indices.size)
+ all_indices = _ramses_reader.get_array_indices_lists(
+ hilbert_indices, unique_indices)
+ for curve_index, my_octs in zip(unique_indices, all_indices):
#print "Handling", curve_index
- my_octs = (hilbert_indices == curve_index)
+ #my_octs = (hilbert_indices == curve_index)
dleft_index = left_index[my_octs,:]
- dright_index = left_index[my_octs,:] + 2
- ddims = (dright_index * 0) + 2
dfl = fl[my_octs,:]
initial_left = na.min(dleft_index, axis=0)
- idims = (na.max(dright_index, axis=0) - initial_left).ravel()
- #if level > 6: insert_ipython()
+ idims = (na.max(dleft_index, axis=0) - initial_left).ravel()+2
+ #if level > 10: insert_ipython()
#print initial_left, idims
psg = _ramses_reader.ProtoSubgrid(initial_left, idims,
- dleft_index, dright_index, ddims, dfl)
+ dleft_index, dfl)
if psg.efficiency <= 0: continue
self.num_deep = 0
psgs.extend(self._recursive_patch_splitting(
psg, idims, initial_left,
- dleft_index, dright_index, ddims, dfl))
+ dleft_index, dfl))
+ print "Done with level % 2i" % (level)
pbar.finish()
self.proto_grids.append(psgs)
sums = na.zeros(3, dtype='int64')
@@ -212,7 +216,7 @@
@num_deep_inc
def _recursive_patch_splitting(self, psg, dims, ind,
- left_index, right_index, gdims, fl):
+ left_index, fl):
min_eff = 0.1 # This isn't always respected.
if self.num_deep > 40:
# If we've recursed more than 100 times, we give up.
@@ -232,13 +236,13 @@
li_l = ind.copy()
if na.any(dims_l <= 0): return [psg]
L = _ramses_reader.ProtoSubgrid(
- li_l, dims_l, left_index, right_index, gdims, fl)
+ li_l, dims_l, left_index, fl)
#print " " * self.num_deep + "L", tt, L.efficiency
if L.efficiency > 1.0: raise RuntimeError
if L.efficiency <= 0.0: L = []
elif L.efficiency < min_eff:
L = self._recursive_patch_splitting(L, dims_l, li_l,
- left_index, right_index, gdims, fl)
+ left_index, fl)
else:
L = [L]
dims_r = dims.copy()
@@ -247,13 +251,13 @@
li_r[ax] += fp
if na.any(dims_r <= 0): return [psg]
R = _ramses_reader.ProtoSubgrid(
- li_r, dims_r, left_index, right_index, gdims, fl)
+ li_r, dims_r, left_index, fl)
#print " " * self.num_deep + "R", tt, R.efficiency
if R.efficiency > 1.0: raise RuntimeError
if R.efficiency <= 0.0: R = []
elif R.efficiency < min_eff:
R = self._recursive_patch_splitting(R, dims_r, li_r,
- left_index, right_index, gdims, fl)
+ left_index, fl)
else:
R = [R]
return L + R
@@ -276,18 +280,16 @@
gi += 1
self.grids = na.array(grids, dtype='object')
- def _get_grid_parents(self, grid, LE, RE):
- mask = na.zeros(self.num_grids, dtype='bool')
- grids, grid_ind = self.get_box_grids(LE, RE)
- mask[grid_ind] = True
- mask = na.logical_and(mask, (self.grid_levels == (grid.Level-1)).flat)
- return self.grids[mask]
-
def _populate_grid_objects(self):
+ mask = na.empty(self.grids.size, dtype='int32')
+ print self.grid_levels.dtype
for gi,g in enumerate(self.grids):
- parents = self._get_grid_parents(g,
- self.grid_left_edge[gi,:],
- self.grid_right_edge[gi,:])
+ get_box_grids_level(self.grid_left_edge[gi,:],
+ self.grid_right_edge[gi,:],
+ g.Level - 1,
+ self.grid_left_edge, self.grid_right_edge,
+ self.grid_levels, mask)
+ parents = self.grids[mask.astype("bool")]
if len(parents) > 0:
g.Parent.extend(parents.tolist())
for p in parents: p.Children.append(g)
--- a/yt/utilities/_amr_utils/misc_utilities.pyx Tue Jun 07 19:03:15 2011 -0400
+++ b/yt/utilities/_amr_utils/misc_utilities.pyx Tue Jun 07 20:22:57 2011 -0400
@@ -56,7 +56,7 @@
int level,
np.ndarray[np.float64_t, ndim=2] left_edges,
np.ndarray[np.float64_t, ndim=2] right_edges,
- np.ndarray[np.int64_t, ndim=2] levels,
+ np.ndarray[np.int32_t, ndim=2] levels,
np.ndarray[np.int32_t, ndim=1] mask):
cdef int i, n
cdef int nx = left_edges.shape[0]
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