[Yt-svn] yt-commit r1147 - trunk/yt/lagos
mturk at wrangler.dreamhost.com
mturk at wrangler.dreamhost.com
Sat Jan 24 15:58:07 PST 2009
Author: mturk
Date: Sat Jan 24 15:58:06 2009
New Revision: 1147
URL: http://yt.spacepope.org/changeset/1147
Log:
Fixed: getting particle vector fields, particle velocity conversion (it was not
being returned in CGS).
Added: ParticleSpinParameter. These agree with enzo_anyl to four decimal
places; I believe the remaining difference to be between using the *actual*
mass inside the virial radius versus 4/3 pi r_vir**3 * rho_vir.
Modified:
trunk/yt/lagos/BaseDataTypes.py
trunk/yt/lagos/DerivedQuantities.py
trunk/yt/lagos/UniversalFields.py
Modified: trunk/yt/lagos/BaseDataTypes.py
==============================================================================
--- trunk/yt/lagos/BaseDataTypes.py (original)
+++ trunk/yt/lagos/BaseDataTypes.py Sat Jan 24 15:58:06 2009
@@ -1277,6 +1277,9 @@
if field in self.pf.field_info and self.pf.field_info[field].particle_type:
if grid.NumberOfParticles == 0: return na.array([])
pointI = self._get_particle_indices(grid)
+ if self.pf.field_info[field].vector_field:
+ f = grid[field]
+ return na.array([f[i,:][pointI] for i in range(3)])
return grid[field][pointI].ravel()
if field in self.pf.field_info and self.pf.field_info[field].vector_field:
pointI = self._get_point_indices(grid)
Modified: trunk/yt/lagos/DerivedQuantities.py
==============================================================================
--- trunk/yt/lagos/DerivedQuantities.py (original)
+++ trunk/yt/lagos/DerivedQuantities.py Sat Jan 24 15:58:06 2009
@@ -186,9 +186,9 @@
This function returns the spin parameter for the baryons, but it uses
the particles in calculating enclosed mass.
"""
+ m_enc = data["CellMassMsun"].sum() + data["ParticleMassMsun"].sum()
am = data["SpecificAngularMomentum"]*data["CellMassMsun"]
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
@@ -201,11 +201,26 @@
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("BaryonSpinParameter", function=_BaryonSpinParameter,
combine_function=_combBaryonSpinParameter, n_ret=4)
+def _ParticleSpinParameter(data):
+ """
+ This function returns the spin parameter for the baryons, but it uses
+ the particles in calculating enclosed mass.
+ """
+ m_enc = data["CellMassMsun"].sum() + data["ParticleMassMsun"].sum()
+ am = data["ParticleSpecificAngularMomentum"]*data["ParticleMassMsun"]
+ if am.size == 0: return (na.zeros((3,), dtype='float64'), m_enc, 0, 0)
+ j_mag = am.sum(axis=1)
+ e_term_pre = na.sum(data["ParticleMassMsun"]
+ *data["ParticleVelocityMagnitude"]**2.0)
+ weight=data["ParticleMassMsun"].sum()
+ return j_mag, m_enc, e_term_pre, weight
+add_quantity("ParticleSpinParameter", function=_ParticleSpinParameter,
+ combine_function=_combBaryonSpinParameter, n_ret=4)
+
def _IsBound(data, truncate = True, include_thermal_energy = False):
"""
This returns whether or not the object is gravitationally bound
Modified: trunk/yt/lagos/UniversalFields.py
==============================================================================
--- trunk/yt/lagos/UniversalFields.py (original)
+++ trunk/yt/lagos/UniversalFields.py Sat Jan 24 15:58:06 2009
@@ -139,12 +139,23 @@
return data._read_data(p_field.replace("_"," ")).astype('float64')
return _Particles
for pf in ["index", "type", "mass"] + \
- ["velocity_%s" % ax for ax in 'xyz'] + \
["position_%s" % ax for ax in 'xyz']:
pfunc = particle_func("particle_%s" % (pf))
add_field("particle_%s" % pf, function=pfunc,
validators = [ValidateSpatial(0)],
particle_type=True)
+def _get_vel_convert(ax):
+ def _convert_p_vel(data):
+ return data.convert("%s-velocity" % ax)
+ return _convert_p_vel
+for ax in 'xyz':
+ pf = "particle_velocity_%s" % ax
+ pfunc = particle_func(pf)
+ cfunc = _get_vel_convert(ax)
+ add_field(pf, function=pfunc, convert_function=cfunc,
+ validators = [ValidateSpatial(0)],
+ particle_type=True)
+
for pf in ["creation_time", "dynamical_time", "metallicity_fraction"]:
pfunc = particle_func(pf)
add_field(pf, function=pfunc,
@@ -203,6 +214,18 @@
convert_function=_convertCourantTimeStep,
units=r"$\rm{s}$")
+def _ParticleVelocityMagnitude(field, data):
+ """M{|v|}"""
+ bulk_velocity = data.get_field_parameter("bulk_velocity")
+ if bulk_velocity == None:
+ bulk_velocity = na.zeros(3)
+ return ( (data["particle_velocity_x"]-bulk_velocity[0])**2.0 + \
+ (data["particle_velocity_y"]-bulk_velocity[1])**2.0 + \
+ (data["particle_velocity_z"]-bulk_velocity[2])**2.0 )**(1.0/2.0)
+add_field("ParticleVelocityMagnitude", function=_ParticleVelocityMagnitude,
+ particle_type=True,
+ take_log=False, units=r"\rm{cm}/\rm{s}")
+
def _VelocityMagnitude(field, data):
"""M{|v|}"""
bulk_velocity = data.get_field_parameter("bulk_velocity")
@@ -530,18 +553,18 @@
v_vec = na.array([xv,yv,zv], dtype='float64')
return na.cross(r_vec, v_vec, axis=0)
add_field("ParticleSpecificAngularMomentum",
- function=_ParticleSpecificAngularMomentum,
+ function=_ParticleSpecificAngularMomentum, particle_type=True,
convert_function=_convertSpecificAngularMomentum, vector_field=True,
units=r"\rm{cm}^2/\rm{s}", validators=[ValidateParameter('center')])
def _convertSpecificAngularMomentumKMSMPC(data):
return data.convert("mpc")/1e5
add_field("ParticleSpecificAngularMomentumKMSMPC",
- function=_ParticleSpecificAngularMomentum,
+ function=_ParticleSpecificAngularMomentum, particle_type=True,
convert_function=_convertSpecificAngularMomentumKMSMPC, vector_field=True,
units=r"\rm{km}\rm{Mpc}/\rm{s}", validators=[ValidateParameter('center')])
def _ParticleAngularMomentum(field, data):
return data["ParticleMass"] * data["ParticleSpecificAngularMomentum"]
-add_field("ParticleAngularMomentum",
+add_field("ParticleAngularMomentum",
function=_ParticleAngularMomentum, units=r"\rm{g}\/\rm{cm}^2/\rm{s}",
particle_type=True)
def _ParticleAngularMomentumMSUNKMSMPC(field, data):
More information about the yt-svn
mailing list