Is anybody here an expert on the way enzo_anyl calculates spin
parameter? I was hoping to convince myself that:
a) I understand how it does it in enzo_anyl
b) That is being replicated in DerivedQuantities.py
I have gone back and forth on this with Brian and Britton, but that
was almost a year ago. I think it needs another look. Here is how yt
works right now:
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()
so we get the sum of Lx, Ly, Lz, the total mass, the total kinetic
energy in the baryons, and the total *baryon* mass.
then during the combination step:
W = weight.sum()
M = m_enc.sum()
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)
What this does it combine all the weights from the individual grid or
processors to get the total baryon mass in the entire region, the
entire *enclosed mass* (which includes the particles), and then the
magnitude of the angular momentum vector for all the enclosed baryons,
which gets divided by the enclosed *baryon mass* to get the average
specific angular momentum for the region. E is then the total kinetic
energy divided by the total enclosed mass, which gives a
characteristic baryon velocity. Finally, we take the average specific
angular momentum, multiply that by the characteristic velocity, and
divide by the total enclosed (baryon+particle) mass.
Does this make sense? Should the characteristic velocity and angular
momentum include the particles? Or does this not matter?
-Matt