[Yt-svn] yt-commit r421 - trunk/yt/lagos
mturk at wrangler.dreamhost.com
mturk at wrangler.dreamhost.com
Thu May 1 08:14:32 PDT 2008
Author: mturk
Date: Thu May 1 08:14:31 2008
New Revision: 421
URL: http://yt.spacepope.org/changeset/421
Log:
Okay, quantities all appear to be working in both lazy and unlazy methods. I
will be comparing against enzo_anyl output shortly.
* Removed dependency on SciPy for finding the maximum location in a
grid/dataset. How did I miss this?
* Particles are now cast to float64 manually. This improved accuracy with the
DerivedQuantities, specifically spin parameter, because of the varying manners
of summing.
Modified:
trunk/yt/lagos/BaseDataTypes.py
trunk/yt/lagos/BaseGridType.py
trunk/yt/lagos/DerivedFields.py
trunk/yt/lagos/DerivedQuantities.py
Modified: trunk/yt/lagos/BaseDataTypes.py
==============================================================================
--- trunk/yt/lagos/BaseDataTypes.py (original)
+++ trunk/yt/lagos/BaseDataTypes.py Thu May 1 08:14:31 2008
@@ -506,7 +506,7 @@
# First we try all three, see which has the best result:
vecs = na.identity(3)
_t = na.cross(self._norm_vec, vecs).sum(axis=1)
- ax = nd.maximum_position(_t)
+ ax = _t.argmax()
self._x_vec = na.cross(vecs[ax,:], self._norm_vec).ravel()
self._x_vec /= na.sqrt(na.dot(self._x_vec, self._x_vec))
self._y_vec = na.cross(self._norm_vec, self._x_vec).ravel()
Modified: trunk/yt/lagos/BaseGridType.py
==============================================================================
--- trunk/yt/lagos/BaseGridType.py (original)
+++ trunk/yt/lagos/BaseGridType.py Thu May 1 08:14:31 2008
@@ -252,7 +252,8 @@
@param field: field to check
@type field: string
"""
- coord=nd.maximum_position(self[field]*self.child_mask)
+ coord1d=(self[field]*self.child_mask).argmax()
+ coord=na.unravel_index(coord1d, self[field].shape)
val = self[field][coord]
return val, coord
@@ -263,8 +264,8 @@
@param field: field to check
@type field: string
"""
- coord=nd.minimum_position(self[field])
- val = self[field][coord]
+ coord1d=(self[field]*self.child_mask).argmin()
+ coord=na.unravel_index(coord1d, self[field].shape)
return val, coord
def get_position(self, coord):
Modified: trunk/yt/lagos/DerivedFields.py
==============================================================================
--- trunk/yt/lagos/DerivedFields.py (original)
+++ trunk/yt/lagos/DerivedFields.py Thu May 1 08:14:31 2008
@@ -285,16 +285,14 @@
def particle_func(p_field):
def _Particles(field, data):
- try:
- particles = data._read_data("particle_%s" % p_field)
- except data._read_exception:
- particles = na.array([], dtype='float64')
- return particles
+ if not data.NumberOfParticles > 0:
+ return na.array([], dtype='float64')
+ return data._read_data(p_field).astype('float64')
return _Particles
for pf in ["index","type"] + \
["velocity_%s" % ax for ax in 'xyz'] + \
["position_%s" % ax for ax in 'xyz']:
- pfunc = particle_func(pf)
+ pfunc = particle_func("particle_%s" % (pf))
add_field("particle_%s" % pf, function=pfunc,
validators = [ValidateSpatial(0)],
particle_type=True)
@@ -302,7 +300,7 @@
validators=[ValidateSpatial(0)], particle_type=True)
def _ParticleMass(field, data):
- particles = data["particle mass"] * \
+ particles = data["particle mass"].astype('float64') * \
just_one(data["CellVolumeCode"].ravel())
# Note that we mandate grid-type here, so this is okay
return particles
@@ -574,10 +572,10 @@
zv = data["z-velocity"]
center = data.get_field_parameter('center')
- coords = na.array([data['x'],data['y'],data['z']])
+ coords = na.array([data['x'],data['y'],data['z']], dtype='float64')
new_shape = tuple([3] + [1]*(len(coords.shape)-1))
r_vec = coords - na.reshape(center,new_shape)
- v_vec = na.array([xv,yv,zv])
+ v_vec = na.array([xv,yv,zv], dtype='float64')
return na.cross(r_vec, v_vec, axis=0)
def _convertSpecificAngularMomentum(data):
return data.convert("cm")
Modified: trunk/yt/lagos/DerivedQuantities.py
==============================================================================
--- trunk/yt/lagos/DerivedQuantities.py (original)
+++ trunk/yt/lagos/DerivedQuantities.py Thu May 1 08:14:31 2008
@@ -50,10 +50,13 @@
def call_func_lazy(args, kwargs):
n_ret = quantities_object.n_ret
retvals = [ [] for i in range(n_ret)]
+ pbar = get_pbar("Calculating ", len(data_object._grids))
for gi,g in enumerate(data_object._grids):
rv = func(GridChildMaskWrapper(g, data_object), *args, **kwargs)
for i in range(n_ret): retvals[i].append(rv[i])
g.clear_data()
+ pbar.update(gi)
+ pbar.finish()
retvals = [na.array(retvals[i]) for i in range(n_ret)]
return c_func(data_object, *retvals)
def call_func(*args, **kwargs):
@@ -114,31 +117,38 @@
add_quantity("WeightedAverageQuantity", function=_WeightedAverageQuantity,
combine_function=_combWeightedAverageQuantity, n_ret = 2)
-"""
-j_vec = na.sum(sp["SpecificAngularMomentum"]*sp["CellMassMsun"],axis=1)/m_weight
-j_mag = na.sqrt(na.sum(j_vec**2.0)) # cm^2 / s
-
-eterm = na.sqrt(0.5*na.sum(sp["CellMassMsun"]*sp["VelocityMagnitude"]**2.0)/m_weight)
-
-G = 6.67e-8 # cm^3 g^-1 s^-2
-m_vir = na.sum(sp["CellMassMsun"]) + na.sum(sp["ParticleMassMsun"])# g
-print "MVIR: %0.5e m_sun" % m_vir
-m_vir *= 1.989e33
-spin = j_mag * eterm / (m_vir * G)
-"""
+def _BulkVelocity(data):
+ xv = (data["x-velocity"] * data["CellMassMsun"]).sum()
+ yv = (data["y-velocity"] * data["CellMassMsun"]).sum()
+ zv = (data["z-velocity"] * data["CellMassMsun"]).sum()
+ w = data["CellMassMsun"].sum()
+ return xv, yv, zv, w
+def _combBulkVelocity(data, xv, yv, zv, w):
+ w = w.sum()
+ xv = xv.sum()/w
+ yv = yv.sum()/w
+ zv = zv.sum()/w
+ return na.array([xv, yv, zv])
+add_quantity("BulkVelocity", function=_BulkVelocity,
+ combine_function=_combBulkVelocity, n_ret=4)
-def _SpinParameter(data):
+def _BaryonSpinParameter(data):
am = data["SpecificAngularMomentum"]*data["CellMassMsun"]
- j_mag = na.sqrt((am**2.0).sum(axis=0)).sum()
- m_enc = data["CellMassMsun"].sum()# + data["ParticleMassMsun"].sum()
- e_term_pre = 0.5*na.sum(data["CellMassMsun"]*data["VelocityMagnitude"]**2.0)
- return j_mag, m_enc, e_term_pre
-def _combSpinParameter(data, j_mag, m_enc, e_term_pre):
+ j_mag = am.sum(axis=1)
+ m_enc = data["CellMassMsun"].sum() + data["ParticleMassMsun"].sum()
+ e_term_pre = na.sum(data["CellMassMsun"]*data["VelocityMagnitude"]**2.0)
+ weight=data["CellMassMsun"].sum()
+ return j_mag, m_enc, e_term_pre, weight
+def _combBaryonSpinParameter(data, j_mag, m_enc, e_term_pre, weight):
+ # Because it's a vector field, we have to ensure we have enough dimensions
+ if len(j_mag.shape) < 2: j_mag = na.expand_dims(j_mag, 0)
+ W = weight.sum()
M = m_enc.sum()
- J = j_mag.sum()/M
- E = na.sqrt(e_term_pre.sum()/M)
+ J = na.sqrt(((j_mag.sum(axis=0))**2.0).sum())/W
+ E = na.sqrt(e_term_pre.sum()/W)
G = 6.67e-8 # cm^3 g^-1 s^-2
spin = J * E / (M*1.989e33*G)
+ print "WEIGHT", W, M, J, E, G, spin
return spin
-add_quantity("SpinParameter", function=_SpinParameter,
- combine_function=_combSpinParameter, n_ret=3)
+add_quantity("BaryonSpinParameter", function=_BaryonSpinParameter,
+ combine_function=_combBaryonSpinParameter, n_ret=4)
More information about the yt-svn
mailing list