List of All Fields

Universal Field List

AbsDivV

  • Units: \rm{s}^{-1}
  • Particle Type: False

Field Source

def _AbsDivV(field, data):
    return na.abs(data['DivV'])

Convert Function Source

No source available.

AngularMomentum

  • Units: \rm{g}\/\rm{cm}^2/\rm{s}
  • Particle Type: False

Field Source

def _AngularMomentum(field, data):
    return data["CellMass"] * data["SpecificAngularMomentum"]

Convert Function Source

No source available.

AngularMomentumMSUNKMSMPC

  • Units: M_{\odot}\rm{km}\rm{Mpc}/\rm{s}
  • Particle Type: False

Field Source

def _AngularMomentum(field, data):
    return data["CellMass"] * data["SpecificAngularMomentum"]

Convert Function Source

No source available.

AngularMomentumX

  • Units: \rm{g}\/\rm{cm}^2/\rm{s}
  • Particle Type: False

Field Source

def _AngularMomentumX(field, data):
    return data["CellMass"] * data["SpecificAngularMomentumX"]

Convert Function Source

No source available.

AngularMomentumY

  • Units: \rm{g}\/\rm{cm}^2/\rm{s}
  • Particle Type: False

Field Source

def _AngularMomentumY(field, data):
    return data["CellMass"] * data["SpecificAngularMomentumY"]

Convert Function Source

No source available.

AngularMomentumZ

  • Units: \rm{g}\/\rm{cm}^2/\rm{s}
  • Particle Type: False

Field Source

def _AngularMomentumZ(field, data):
    return data["CellMass"] * data["SpecificAngularMomentumZ"]

Convert Function Source

No source available.

AveragedDensity

  • Particle Type: False

Field Source

def _AveragedDensity(field, data):
    nx, ny, nz = data["Density"].shape
    new_field = na.zeros((nx-2,ny-2,nz-2), dtype='float64')
    weight_field = na.zeros((nx-2,ny-2,nz-2), dtype='float64')
    i_i, j_i, k_i = na.mgrid[0:3,0:3,0:3]
    for i,j,k in zip(i_i.ravel(),j_i.ravel(),k_i.ravel()):
        sl = [slice(i,nx-(2-i)),slice(j,ny-(2-j)),slice(k,nz-(2-k))]
        new_field += data["Density"][sl] * data["CellMass"][sl]
        weight_field += data["CellMass"][sl]
    # Now some fancy footwork
    new_field2 = na.zeros((nx,ny,nz))
    new_field2[1:-1,1:-1,1:-1] = new_field/weight_field
    return new_field2

Convert Function Source

No source available.

CellMass

  • Units: \rm{g}
  • Particle Type: False

Field Source

def _CellMass(field, data):
    return data["Density"] * data["CellVolume"]

Convert Function Source

No source available.

CellMassCode

  • Particle Type: False

Field Source

def _CellMassCode(field, data):
    return data["Density"] * data["CellVolumeCode"]

Convert Function Source

def _convertCellMassCode(data):
    return 1.0/data.convert("Density")

CellMassMsun

  • Units: M_{\odot}
  • Particle Type: False

Field Source

def _CellMass(field, data):
    return data["Density"] * data["CellVolume"]

Convert Function Source

def _convertCellMassMsun(data):
    return 5.027854e-34 # g^-1

CellVolume

  • Units: \rm{cm}^3
  • Particle Type: False

Field Source

def _CellVolume(field, data):
    if data['dx'].size == 1:
        try:
            return data['dx']*data['dy']*data['dx']*\
                na.ones(data.ActiveDimensions, dtype='float64')
        except AttributeError:
            return data['dx']*data['dy']*data['dx']
    return data["dx"]*data["dy"]*data["dz"]

Convert Function Source

def _ConvertCellVolumeCGS(data):
    return data.convert("cm")**3.0

CellVolumeCode

  • Units: \rm{BoxVolume}^3
  • Particle Type: False

Field Source

def _CellVolume(field, data):
    if data['dx'].size == 1:
        try:
            return data['dx']*data['dy']*data['dx']*\
                na.ones(data.ActiveDimensions, dtype='float64')
        except AttributeError:
            return data['dx']*data['dy']*data['dx']
    return data["dx"]*data["dy"]*data["dz"]

Convert Function Source

No source available.

CellVolumeMpc

  • Units: \rm{Mpc}^3
  • Particle Type: False

Field Source

def _CellVolume(field, data):
    if data['dx'].size == 1:
        try:
            return data['dx']*data['dy']*data['dx']*\
                na.ones(data.ActiveDimensions, dtype='float64')
        except AttributeError:
            return data['dx']*data['dy']*data['dx']
    return data["dx"]*data["dy"]*data["dz"]

Convert Function Source

def _ConvertCellVolumeMpc(data):
    return data.convert("mpc")**3.0

CellsPerBin

  • Particle Type: False

Field Source

def _Ones(field, data):
    return na.ones(data.ActiveDimensions, dtype='float64')

Convert Function Source

No source available.

Contours

  • Particle Type: False

Field Source

def _Contours(field, data):
    return na.ones(data["Density"].shape)*-1

Convert Function Source

No source available.

CourantTimeStep

  • Units: $\rm{s}$
  • Particle Type: False

Field Source

def _CourantTimeStep(field, data):
    t1 = data['dx'] / (
        data["SoundSpeed"] + \
        abs(data["x-velocity"]))
    t2 = data['dy'] / (
        data["SoundSpeed"] + \
        abs(data["y-velocity"]))
    t3 = data['dz'] / (
        data["SoundSpeed"] + \
        abs(data["z-velocity"]))
    return na.minimum(na.minimum(t1,t2),t3)

Convert Function Source

def _convertCourantTimeStep(data):
    # SoundSpeed and z-velocity are in cm/s, dx is in code
    return data.convert("cm")

CuttingPlaneVelocityX

  • Units: \rm{km}/\rm{s}
  • Particle Type: False

Field Source

def _CuttingPlaneVelocityX(field, data):
    x_vec, y_vec, z_vec = [data.get_field_parameter("cp_%s_vec" % (ax))
                           for ax in 'xyz']
    bulk_velocity = data.get_field_parameter("bulk_velocity")
    if bulk_velocity == None:
        bulk_velocity = na.zeros(3)
    v_vec = na.array([data["%s-velocity" % ax] for ax in 'xyz']) \
                - bulk_velocity[...,na.newaxis]
    return na.dot(x_vec, v_vec)

Convert Function Source

No source available.

CuttingPlaneVelocityY

  • Units: \rm{km}/\rm{s}
  • Particle Type: False

Field Source

def _CuttingPlaneVelocityY(field, data):
    x_vec, y_vec, z_vec = [data.get_field_parameter("cp_%s_vec" % (ax))
                           for ax in 'xyz']
    bulk_velocity = data.get_field_parameter("bulk_velocity")
    if bulk_velocity == None:
        bulk_velocity = na.zeros(3)
    v_vec = na.array([data["%s-velocity" % ax] for ax in 'xyz']) \
                - bulk_velocity[...,na.newaxis]
    return na.dot(y_vec, v_vec)

Convert Function Source

No source available.

DiskAngle

  • Particle Type: False

Field Source

def _DiskAngle(field, data):
    # We make both r_vec and h_vec into unit vectors
    center = data.get_field_parameter("center")
    r_vec = na.array([data["x"] - center[0],
                      data["y"] - center[1],
                      data["z"] - center[2]])
    r_vec = r_vec/na.sqrt((r_vec**2.0).sum(axis=0))
    h_vec = na.array(data.get_field_parameter("height_vector"))
    dp = r_vec[0,:] * h_vec[0] \
       + r_vec[1,:] * h_vec[1] \
       + r_vec[2,:] * h_vec[2]
    return na.arccos(dp)

Convert Function Source

No source available.

DivV

  • Units: \rm{s}^{-1}
  • Particle Type: False

Field Source

def _DivV(field, data):
    # We need to set up stencils
    if data.pf["HydroMethod"] == 2:
        sl_left = slice(None,-2,None)
        sl_right = slice(1,-1,None)
        div_fac = 1.0
    else:
        sl_left = slice(None,-2,None)
        sl_right = slice(2,None,None)
        div_fac = 2.0
    ds = div_fac * data['dx'].flat[0]
    f  = data["x-velocity"][sl_right,1:-1,1:-1]/ds
    f -= data["x-velocity"][sl_left ,1:-1,1:-1]/ds
    ds = div_fac * data['dy'].flat[0]
    f += data["y-velocity"][1:-1,sl_right,1:-1]/ds
    f -= data["y-velocity"][1:-1,sl_left ,1:-1]/ds
    ds = div_fac * data['dz'].flat[0]
    f += data["z-velocity"][1:-1,1:-1,sl_right]/ds
    f -= data["z-velocity"][1:-1,1:-1,sl_left ]/ds
    new_field = na.zeros(data["x-velocity"].shape, dtype='float64')
    new_field[1:-1,1:-1,1:-1] = f
    return new_field

Convert Function Source

def _convertDivV(data):
    return data.convert("cm")**-1.0

DynamicalTime

  • Units: \rm{s}
  • Particle Type: False

Field Source

def _DynamicalTime(field, data):
    """
    The formulation for the dynamical time is:
    M{sqrt(3pi/(16*G*rho))} or M{sqrt(3pi/(16G))*rho^-(1/2)}
    Note that we return in our natural units already
    """
    return (3.0*pi/(16*G*data["Density"]))**(1./2.)

Convert Function Source

No source available.

Entropy

  • Units: \rm{ergs}\/\rm{cm}^{2}
  • Particle Type: False

Field Source

def _Entropy(field, data):
    return (kboltz/mh) * data["Temperature"] / \
           (data["MeanMolecularWeight"] * data["Density"]**(2./3.))

Convert Function Source

No source available.

GridIndices

  • Particle Type: False

Field Source

def _GridIndices(field, data):
    return na.ones(data["Ones"].shape)*(data.id-data._id_offset)

Convert Function Source

No source available.

GridLevel

  • Particle Type: False

Field Source

def _GridLevel(field, data):
    return na.ones(data.ActiveDimensions)*(data.Level)

Convert Function Source

No source available.

Height

  • Units: cm
  • Particle Type: False

Field Source

def _Height(field, data):
    # We take the dot product of the radius vector with the height-vector
    center = data.get_field_parameter("center")
    r_vec = na.array([data["x"] - center[0],
                      data["y"] - center[1],
                      data["z"] - center[2]])
    h_vec = na.array(data.get_field_parameter("height_vector"))
    h_vec = h_vec / na.sqrt(h_vec[0]**2.0+
                            h_vec[1]**2.0+
                            h_vec[2]**2.0)
    height = r_vec[0,:] * h_vec[0] \
           + r_vec[1,:] * h_vec[1] \
           + r_vec[2,:] * h_vec[2]
    return na.abs(height)

Convert Function Source

def _convertHeight(data):
    return data.convert("cm")

HeightAU

  • Units: AU
  • Particle Type: False

Field Source

def _Height(field, data):
    # We take the dot product of the radius vector with the height-vector
    center = data.get_field_parameter("center")
    r_vec = na.array([data["x"] - center[0],
                      data["y"] - center[1],
                      data["z"] - center[2]])
    h_vec = na.array(data.get_field_parameter("height_vector"))
    h_vec = h_vec / na.sqrt(h_vec[0]**2.0+
                            h_vec[1]**2.0+
                            h_vec[2]**2.0)
    height = r_vec[0,:] * h_vec[0] \
           + r_vec[1,:] * h_vec[1] \
           + r_vec[2,:] * h_vec[2]
    return na.abs(height)

Convert Function Source

def _convertHeightAU(data):
    return data.convert("au")

JeansMassMsun

  • Units: \rm{M_{\odot}}
  • Particle Type: False

Field Source

def _JeansMassMsun(field,data):
    MJ_constant = (((5*kboltz)/(G*mh))**(1.5)) * \
    (3/(4*3.1415926535897931))**(0.5) / 1.989e33

    return (MJ_constant *
            ((data["Temperature"]/data["MeanMolecularWeight"])**(1.5)) *
            (data["Density"]**(-0.5)))

Convert Function Source

No source available.

MachNumber

  • Particle Type: False

Field Source

def _MachNumber(field, data):
    """M{|v|/t_sound}"""
    return data["VelocityMagnitude"] / data["SoundSpeed"]

Convert Function Source

No source available.

Matter_Density

  • Units: \rm{g}/\rm{cm^3}
  • Particle Type: False

Field Source

def _Matter_Density(field,data):
    return (data['Density'] + data['Dark_Matter_Density'])

Convert Function Source

No source available.

MeanMolecularWeight

  • Particle Type: False

Field Source

def _MeanMolecularWeight(field,data):
    return (data["Density"] / (mh *data["NumberDensity"]))

Convert Function Source

No source available.

Ones

  • Particle Type: False

Field Source

def _Ones(field, data):
    return na.ones(data.ActiveDimensions, dtype='float64')

Convert Function Source

No source available.

OnesOverDx

  • Particle Type: False

Field Source

def _OnesOverDx(field, data):
    return na.ones(data["Ones"].shape,
                   dtype=data["Density"].dtype)/data['dx']

Convert Function Source

No source available.

ParticleAge

  • Particle Type: True

Field Source

def _ParticleAge(field, data):
    current_time = data.pf["InitialTime"]
    return (current_time - data["creation_time"])

Convert Function Source

def _convertParticleAge(data):
    return data.convert("years")

ParticleAngularMomentum

  • Units: \rm{g}\/\rm{cm}^2/\rm{s}
  • Particle Type: True

Field Source

def _ParticleAngularMomentum(field, data):
    return data["ParticleMass"] * data["ParticleSpecificAngularMomentum"]

Convert Function Source

No source available.

ParticleAngularMomentumMSUNKMSMPC

  • Units: M_{\odot}\rm{km}\rm{Mpc}/\rm{s}
  • Particle Type: True

Field Source

def _ParticleAngularMomentumMSUNKMSMPC(field, data):
    return data["ParticleMass"] * data["ParticleSpecificAngularMomentumKMSMPC"]

Convert Function Source

No source available.

ParticleMass

  • Particle Type: True

Field Source

def _ParticleMass(field, data):
    particles = data["particle_mass"].astype('float64') * \
                just_one(data["CellVolumeCode"].ravel())
    # Note that we mandate grid-type here, so this is okay
    return particles

Convert Function Source

def _convertParticleMass(data):
    return data.convert("Density")*(data.convert("cm")**3.0)

ParticleMassMsun

  • Particle Type: True

Field Source

def _ParticleMass(field, data):
    particles = data["particle_mass"].astype('float64') * \
                just_one(data["CellVolumeCode"].ravel())
    # Note that we mandate grid-type here, so this is okay
    return particles

Convert Function Source

def _convertParticleMassMsun(data):
    return data.convert("Density")*((data.convert("cm")**3.0)/1.989e33)

ParticleRadius

  • Units: \rm{cm}
  • Particle Type: True

Field Source

def _ParticleRadius(field, data):
    center = data.get_field_parameter("center")
    DW = data.pf["DomainRightEdge"] - data.pf["DomainLeftEdge"]
    radius = na.zeros(data["particle_position_x"].shape, dtype='float64')
    for i, ax in enumerate('xyz'):
        r = na.abs(data["particle_position_%s" % ax] - center[i])
        radius += na.minimum(r, na.abs(DW[i]-r))**2.0
    na.sqrt(radius, radius)
    return radius

Convert Function Source

def _ConvertRadiusCGS(data):
    return data.convert("cm")

ParticleRadiusAU

  • Units: \rm{AU}
  • Particle Type: True

Field Source

def _ParticleRadius(field, data):
    center = data.get_field_parameter("center")
    DW = data.pf["DomainRightEdge"] - data.pf["DomainLeftEdge"]
    radius = na.zeros(data["particle_position_x"].shape, dtype='float64')
    for i, ax in enumerate('xyz'):
        r = na.abs(data["particle_position_%s" % ax] - center[i])
        radius += na.minimum(r, na.abs(DW[i]-r))**2.0
    na.sqrt(radius, radius)
    return radius

Convert Function Source

def _ConvertRadiusAU(data):
    return data.convert("au")

ParticleRadiusCode

  • Particle Type: True

Field Source

def _ParticleRadius(field, data):
    center = data.get_field_parameter("center")
    DW = data.pf["DomainRightEdge"] - data.pf["DomainLeftEdge"]
    radius = na.zeros(data["particle_position_x"].shape, dtype='float64')
    for i, ax in enumerate('xyz'):
        r = na.abs(data["particle_position_%s" % ax] - center[i])
        radius += na.minimum(r, na.abs(DW[i]-r))**2.0
    na.sqrt(radius, radius)
    return radius

Convert Function Source

No source available.

ParticleRadiusMpc

  • Units: \rm{Mpc}
  • Particle Type: True

Field Source

def _ParticleRadius(field, data):
    center = data.get_field_parameter("center")
    DW = data.pf["DomainRightEdge"] - data.pf["DomainLeftEdge"]
    radius = na.zeros(data["particle_position_x"].shape, dtype='float64')
    for i, ax in enumerate('xyz'):
        r = na.abs(data["particle_position_%s" % ax] - center[i])
        radius += na.minimum(r, na.abs(DW[i]-r))**2.0
    na.sqrt(radius, radius)
    return radius

Convert Function Source

def _ConvertRadiusMpc(data):
    return data.convert("mpc")

ParticleRadiuskpc

  • Units: \rm{kpc}
  • Particle Type: True

Field Source

def _ParticleRadius(field, data):
    center = data.get_field_parameter("center")
    DW = data.pf["DomainRightEdge"] - data.pf["DomainLeftEdge"]
    radius = na.zeros(data["particle_position_x"].shape, dtype='float64')
    for i, ax in enumerate('xyz'):
        r = na.abs(data["particle_position_%s" % ax] - center[i])
        radius += na.minimum(r, na.abs(DW[i]-r))**2.0
    na.sqrt(radius, radius)
    return radius

Convert Function Source

def _ConvertRadiuskpc(data):
    return data.convert("kpc")

ParticleRadiuskpch

  • Units: \rm{kpc}/\rm{h}
  • Particle Type: True

Field Source

def _ParticleRadius(field, data):
    center = data.get_field_parameter("center")
    DW = data.pf["DomainRightEdge"] - data.pf["DomainLeftEdge"]
    radius = na.zeros(data["particle_position_x"].shape, dtype='float64')
    for i, ax in enumerate('xyz'):
        r = na.abs(data["particle_position_%s" % ax] - center[i])
        radius += na.minimum(r, na.abs(DW[i]-r))**2.0
    na.sqrt(radius, radius)
    return radius

Convert Function Source

def _ConvertRadiuskpc(data):
    return data.convert("kpc")

ParticleRadiuspc

  • Units: \rm{pc}
  • Particle Type: True

Field Source

def _ParticleRadius(field, data):
    center = data.get_field_parameter("center")
    DW = data.pf["DomainRightEdge"] - data.pf["DomainLeftEdge"]
    radius = na.zeros(data["particle_position_x"].shape, dtype='float64')
    for i, ax in enumerate('xyz'):
        r = na.abs(data["particle_position_%s" % ax] - center[i])
        radius += na.minimum(r, na.abs(DW[i]-r))**2.0
    na.sqrt(radius, radius)
    return radius

Convert Function Source

def _ConvertRadiuspc(data):
    return data.convert("pc")

ParticleSpecificAngularMomentum

  • Units: \rm{cm}^2/\rm{s}
  • Particle Type: True

Field Source

def _ParticleSpecificAngularMomentum(field, data):
    """
    Calculate the angular of a particle velocity.  Returns a vector for each
    particle.
    """
    if data.has_field_parameter("bulk_velocity"):
        bv = data.get_field_parameter("bulk_velocity")
    else: bv = na.zeros(3, dtype='float64')
    xv = data["particle_velocity_x"] - bv[0]
    yv = data["particle_velocity_y"] - bv[1]
    zv = data["particle_velocity_z"] - bv[2]
    center = data.get_field_parameter('center')
    coords = na.array([data['particle_position_x'],
                       data['particle_position_y'],
                       data['particle_position_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], dtype='float64')
    return na.cross(r_vec, v_vec, axis=0)

Convert Function Source

def _convertSpecificAngularMomentum(data):
    return data.convert("cm")

ParticleSpecificAngularMomentumKMSMPC

  • Units: \rm{km}\rm{Mpc}/\rm{s}
  • Particle Type: True

Field Source

def _ParticleSpecificAngularMomentum(field, data):
    """
    Calculate the angular of a particle velocity.  Returns a vector for each
    particle.
    """
    if data.has_field_parameter("bulk_velocity"):
        bv = data.get_field_parameter("bulk_velocity")
    else: bv = na.zeros(3, dtype='float64')
    xv = data["particle_velocity_x"] - bv[0]
    yv = data["particle_velocity_y"] - bv[1]
    zv = data["particle_velocity_z"] - bv[2]
    center = data.get_field_parameter('center')
    coords = na.array([data['particle_position_x'],
                       data['particle_position_y'],
                       data['particle_position_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], dtype='float64')
    return na.cross(r_vec, v_vec, axis=0)

Convert Function Source

def _convertSpecificAngularMomentumKMSMPC(data):
    return data.convert("mpc")/1e5

ParticleVelocityMagnitude

  • Units: \rm{cm}/\rm{s}
  • Particle Type: True

Field Source

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)

Convert Function Source

No source available.

Pressure

  • Units: \rm{dyne}/\rm{cm}^{2}
  • Particle Type: False

Field Source

def _Pressure(field, data):
    """M{(Gamma-1.0)*rho*E}"""
    return (data.pf["Gamma"] - 1.0) * \
           data["Density"] * data["ThermalEnergy"]

Convert Function Source

No source available.

RadialMachNumber

  • Particle Type: False

Field Source

def _RadialMachNumber(field, data):
    """M{|v|/t_sound}"""
    return na.abs(data["RadialVelocity"]) / data["SoundSpeed"]

Convert Function Source

No source available.

RadialVelocity

  • Units: \rm{cm}/\rm{s}
  • Particle Type: False

Field Source

def _RadialVelocity(field, data):
    center = data.get_field_parameter("center")
    bulk_velocity = data.get_field_parameter("bulk_velocity")
    if bulk_velocity == None:
        bulk_velocity = na.zeros(3)
    new_field = ( (data['x']-center[0])*(data["x-velocity"]-bulk_velocity[0])
                + (data['y']-center[1])*(data["y-velocity"]-bulk_velocity[1])
                + (data['z']-center[2])*(data["z-velocity"]-bulk_velocity[2])
                )/data["RadiusCode"]
    if na.any(na.isnan(new_field)): # to fix center = point
        new_field[na.isnan(new_field)] = 0.0
    return new_field

Convert Function Source

No source available.

RadialVelocityABS

  • Units: \rm{cm}/\rm{s}
  • Particle Type: False

Field Source

def _RadialVelocityABS(field, data):
    return na.abs(_RadialVelocity(field, data))

Convert Function Source

No source available.

RadialVelocityKMS

  • Units: \rm{km}/\rm{s}
  • Particle Type: False

Field Source

def _RadialVelocity(field, data):
    center = data.get_field_parameter("center")
    bulk_velocity = data.get_field_parameter("bulk_velocity")
    if bulk_velocity == None:
        bulk_velocity = na.zeros(3)
    new_field = ( (data['x']-center[0])*(data["x-velocity"]-bulk_velocity[0])
                + (data['y']-center[1])*(data["y-velocity"]-bulk_velocity[1])
                + (data['z']-center[2])*(data["z-velocity"]-bulk_velocity[2])
                )/data["RadiusCode"]
    if na.any(na.isnan(new_field)): # to fix center = point
        new_field[na.isnan(new_field)] = 0.0
    return new_field

Convert Function Source

def _ConvertRadialVelocityKMS(data):
    return 1e-5

RadialVelocityKMSABS

  • Units: \rm{km}/\rm{s}
  • Particle Type: False

Field Source

def _RadialVelocityABS(field, data):
    return na.abs(_RadialVelocity(field, data))

Convert Function Source

def _ConvertRadialVelocityKMS(data):
    return 1e-5

Radius

  • Units: \rm{cm}
  • Particle Type: False

Field Source

def _Radius(field, data):
    center = data.get_field_parameter("center")
    DW = data.pf["DomainRightEdge"] - data.pf["DomainLeftEdge"]
    radius = na.zeros(data["x"].shape, dtype='float64')
    for i, ax in enumerate('xyz'):
        r = na.abs(data[ax] - center[i])
        radius += na.minimum(r, na.abs(DW[i]-r))**2.0
    na.sqrt(radius, radius)
    return radius

Convert Function Source

def _ConvertRadiusCGS(data):
    return data.convert("cm")

RadiusAU

  • Units: \rm{AU}
  • Particle Type: False

Field Source

def _Radius(field, data):
    center = data.get_field_parameter("center")
    DW = data.pf["DomainRightEdge"] - data.pf["DomainLeftEdge"]
    radius = na.zeros(data["x"].shape, dtype='float64')
    for i, ax in enumerate('xyz'):
        r = na.abs(data[ax] - center[i])
        radius += na.minimum(r, na.abs(DW[i]-r))**2.0
    na.sqrt(radius, radius)
    return radius

Convert Function Source

def _ConvertRadiusAU(data):
    return data.convert("au")

RadiusCode

  • Particle Type: False

Field Source

def _Radius(field, data):
    center = data.get_field_parameter("center")
    DW = data.pf["DomainRightEdge"] - data.pf["DomainLeftEdge"]
    radius = na.zeros(data["x"].shape, dtype='float64')
    for i, ax in enumerate('xyz'):
        r = na.abs(data[ax] - center[i])
        radius += na.minimum(r, na.abs(DW[i]-r))**2.0
    na.sqrt(radius, radius)
    return radius

Convert Function Source

No source available.

RadiusMpc

  • Units: \rm{Mpc}
  • Particle Type: False

Field Source

def _Radius(field, data):
    center = data.get_field_parameter("center")
    DW = data.pf["DomainRightEdge"] - data.pf["DomainLeftEdge"]
    radius = na.zeros(data["x"].shape, dtype='float64')
    for i, ax in enumerate('xyz'):
        r = na.abs(data[ax] - center[i])
        radius += na.minimum(r, na.abs(DW[i]-r))**2.0
    na.sqrt(radius, radius)
    return radius

Convert Function Source

def _ConvertRadiusMpc(data):
    return data.convert("mpc")

Radiuskpc

  • Units: \rm{kpc}
  • Particle Type: False

Field Source

def _Radius(field, data):
    center = data.get_field_parameter("center")
    DW = data.pf["DomainRightEdge"] - data.pf["DomainLeftEdge"]
    radius = na.zeros(data["x"].shape, dtype='float64')
    for i, ax in enumerate('xyz'):
        r = na.abs(data[ax] - center[i])
        radius += na.minimum(r, na.abs(DW[i]-r))**2.0
    na.sqrt(radius, radius)
    return radius

Convert Function Source

def _ConvertRadiuskpc(data):
    return data.convert("kpc")

Radiuskpch

  • Units: \rm{kpc}/\rm{h}
  • Particle Type: False

Field Source

def _Radius(field, data):
    center = data.get_field_parameter("center")
    DW = data.pf["DomainRightEdge"] - data.pf["DomainLeftEdge"]
    radius = na.zeros(data["x"].shape, dtype='float64')
    for i, ax in enumerate('xyz'):
        r = na.abs(data[ax] - center[i])
        radius += na.minimum(r, na.abs(DW[i]-r))**2.0
    na.sqrt(radius, radius)
    return radius

Convert Function Source

def _ConvertRadiuskpc(data):
    return data.convert("kpc")

Radiuspc

  • Units: \rm{pc}
  • Particle Type: False

Field Source

def _Radius(field, data):
    center = data.get_field_parameter("center")
    DW = data.pf["DomainRightEdge"] - data.pf["DomainLeftEdge"]
    radius = na.zeros(data["x"].shape, dtype='float64')
    for i, ax in enumerate('xyz'):
        r = na.abs(data[ax] - center[i])
        radius += na.minimum(r, na.abs(DW[i]-r))**2.0
    na.sqrt(radius, radius)
    return radius

Convert Function Source

def _ConvertRadiuspc(data):
    return data.convert("pc")

SZKinetic

  • Particle Type: False

Field Source

def _SZKinetic(field, data):
    vel_axis = data.get_field_parameter('axis')
    if vel_axis > 2:
        raise NeedsParameter(['axis'])
    vel = data["%s-velocity" % ({0:'x',1:'y',2:'z'}[vel_axis])]
    return (vel*data["Density"])

Convert Function Source

def _convertSZKinetic(data):
    return 0.88*((sigma_thompson/mh)/clight)

SZY

  • Particle Type: False

Field Source

def _SZY(field, data):
    return (data["Density"]*data["Temperature"])

Convert Function Source

def _convertSZY(data):
    conv = (0.88/mh) * (kboltz)/(me * clight*clight) * sigma_thompson
    return conv

SoundSpeed

  • Units: \rm{cm}/\rm{s}
  • Particle Type: False

Field Source

def _SoundSpeed(field, data):
    if data.pf["EOSType"] == 1:
        return na.ones(data["Density"].shape, dtype='float64') * \
                data.pf["EOSSoundSpeed"]
    return ( data.pf["Gamma"]*data["Pressure"] / \
             data["Density"] )**(1.0/2.0)

Convert Function Source

No source available.

SpecificAngularMomentum

  • Units: \rm{cm}^2/\rm{s}
  • Particle Type: False

Field Source

def _SpecificAngularMomentum(field, data):
    """
    Calculate the angular velocity.  Returns a vector for each cell.
    """
    r_vec = obtain_rvec(data)
    xv, yv, zv = obtain_velocities(data)
    v_vec = na.array([xv,yv,zv], dtype='float64')
    return na.cross(r_vec, v_vec, axis=0)

Convert Function Source

def _convertSpecificAngularMomentum(data):
    return data.convert("cm")

SpecificAngularMomentumKMSMPC

  • Units: \rm{km}\rm{Mpc}/\rm{s}
  • Particle Type: False

Field Source

def _SpecificAngularMomentum(field, data):
    """
    Calculate the angular velocity.  Returns a vector for each cell.
    """
    r_vec = obtain_rvec(data)
    xv, yv, zv = obtain_velocities(data)
    v_vec = na.array([xv,yv,zv], dtype='float64')
    return na.cross(r_vec, v_vec, axis=0)

Convert Function Source

def _convertSpecificAngularMomentumKMSMPC(data):
    return data.convert("mpc")/1e5

SpecificAngularMomentumX

  • Units: \rm{cm}^2/\rm{s}
  • Particle Type: False

Field Source

def _SpecificAngularMomentumX(field, data):
    xv, yv, zv = obtain_velocities(data)
    rv = obtain_rvec(data)
    return yv*rv[2,:] - zv*rv[1,:]

Convert Function Source

def _convertSpecificAngularMomentum(data):
    return data.convert("cm")

SpecificAngularMomentumY

  • Units: \rm{cm}^2/\rm{s}
  • Particle Type: False

Field Source

def _SpecificAngularMomentumY(field, data):
    xv, yv, zv = obtain_velocities(data)
    rv = obtain_rvec(data)
    return -(xv*rv[2,:] - zv*rv[0,:])

Convert Function Source

def _convertSpecificAngularMomentum(data):
    return data.convert("cm")

SpecificAngularMomentumZ

  • Units: \rm{cm}^2/\rm{s}
  • Particle Type: False

Field Source

def _SpecificAngularMomentumZ(field, data):
    xv, yv, zv = obtain_velocities(data)
    rv = obtain_rvec(data)
    return xv*rv[1,:] - yv*rv[0,:]

Convert Function Source

def _convertSpecificAngularMomentum(data):
    return data.convert("cm")

StarMassMsun

  • Units: M_{\odot}
  • Particle Type: False

Field Source

def _StarMass(field,data):
    return data["star_density_pyx"] * data["CellVolume"]

Convert Function Source

def _convertCellMassMsun(data):
    return 5.027854e-34 # g^-1

TangentialOverVelocityMagnitude

  • Particle Type: False

Field Source

def _TangentialOverVelocityMagnitude(field, data):
    return na.abs(data["TangentialVelocity"])/na.abs(data["VelocityMagnitude"])

Convert Function Source

No source available.

TangentialVelocity

  • Units: \rm{cm}/\rm{s}
  • Particle Type: False

Field Source

def _TangentialVelocity(field, data):
    return na.sqrt(data["VelocityMagnitude"]**2.0
                 - data["RadialVelocity"]**2.0)

Convert Function Source

No source available.

TotalMassMsun

  • Units: M_{\odot}
  • Particle Type: False

Field Source

def _TotalMass(field,data):
    return (data["Density"]+data["Dark_Matter_Density"]) * data["CellVolume"]

Convert Function Source

def _convertCellMassMsun(data):
    return 5.027854e-34 # g^-1

VelocityMagnitude

  • Units: \rm{cm}/\rm{s}
  • Particle Type: False

Field Source

def _VelocityMagnitude(field, data):
    """M{|v|}"""
    bulk_velocity = data.get_field_parameter("bulk_velocity")
    if bulk_velocity == None:
        bulk_velocity = na.zeros(3)
    return ( (data["x-velocity"]-bulk_velocity[0])**2.0 + \
             (data["y-velocity"]-bulk_velocity[1])**2.0 + \
             (data["z-velocity"]-bulk_velocity[2])**2.0 )**(1.0/2.0)

Convert Function Source

No source available.

XRayEmissivity

  • Particle Type: False

Field Source

def _XRayEmissivity(field, data):
    return ((data["Density"].astype('float64')**2.0) \
            *data["Temperature"]**0.5)

Convert Function Source

def _convertXRayEmissivity(data):
    return 2.168e60

creation_time

  • Particle Type: True

Field Source

def _Particles(field, data):
    io = data.hierarchy.io
    if not data.NumberOfParticles > 0:
        return na.array([], dtype=dtype)
    try:
        return io._read_data_set(data, p_field).astype(dtype)
    except io._read_exception:
        pass
    # This is bad.  But it's the best idea I have right now.
    return data._read_data(p_field.replace("_"," ")).astype(dtype)

Convert Function Source

No source available.

dx

  • Particle Type: False

Field Source

def _dx(field, data):
    return data.dds[0]
    return na.ones(data.ActiveDimensions, dtype='float64') * data.dds[0]

Convert Function Source

No source available.

dy

  • Particle Type: False

Field Source

def _dy(field, data):
    return data.dds[1]
    return na.ones(data.ActiveDimensions, dtype='float64') * data.dds[1]

Convert Function Source

No source available.

dynamical_time

  • Particle Type: True

Field Source

def _Particles(field, data):
    io = data.hierarchy.io
    if not data.NumberOfParticles > 0:
        return na.array([], dtype=dtype)
    try:
        return io._read_data_set(data, p_field).astype(dtype)
    except io._read_exception:
        pass
    # This is bad.  But it's the best idea I have right now.
    return data._read_data(p_field.replace("_"," ")).astype(dtype)

Convert Function Source

No source available.

dz

  • Particle Type: False

Field Source

def _dz(field, data):
    return data.dds[2]
    return na.ones(data.ActiveDimensions, dtype='float64') * data.dds[2]

Convert Function Source

No source available.

metallicity_fraction

  • Particle Type: True

Field Source

def _Particles(field, data):
    io = data.hierarchy.io
    if not data.NumberOfParticles > 0:
        return na.array([], dtype=dtype)
    try:
        return io._read_data_set(data, p_field).astype(dtype)
    except io._read_exception:
        pass
    # This is bad.  But it's the best idea I have right now.
    return data._read_data(p_field.replace("_"," ")).astype(dtype)

Convert Function Source

No source available.

particle_density_pyx

  • Particle Type: False

Field Source

def _pdensity_pyx(field, data):
    blank = na.zeros(data.ActiveDimensions, dtype='float32')
    if data.NumberOfParticles == 0: return blank
    CICDeposit_3(data["particle_position_x"].astype(na.float64),
                 data["particle_position_y"].astype(na.float64),
                 data["particle_position_z"].astype(na.float64),
                 data["particle_mass"].astype(na.float32),
                 na.int64(data.NumberOfParticles),
                 blank, na.array(data.LeftEdge).astype(na.float64),
                 na.array(data.ActiveDimensions).astype(na.int32),
                 na.float64(data['dx']))
    return blank

Convert Function Source

def _convertDensity(data):
    return data.convert("Density")

particle_index

  • Particle Type: True

Field Source

def _Particles(field, data):
    io = data.hierarchy.io
    if not data.NumberOfParticles > 0:
        return na.array([], dtype=dtype)
    try:
        return io._read_data_set(data, p_field).astype(dtype)
    except io._read_exception:
        pass
    # This is bad.  But it's the best idea I have right now.
    return data._read_data(p_field.replace("_"," ")).astype(dtype)

Convert Function Source

def _convRetainInt(data):
    return 1

particle_mass

  • Particle Type: True

Field Source

def _Particles(field, data):
    io = data.hierarchy.io
    if not data.NumberOfParticles > 0:
        return na.array([], dtype=dtype)
    try:
        return io._read_data_set(data, p_field).astype(dtype)
    except io._read_exception:
        pass
    # This is bad.  But it's the best idea I have right now.
    return data._read_data(p_field.replace("_"," ")).astype(dtype)

Convert Function Source

No source available.

particle_position_x

  • Particle Type: True

Field Source

def _Particles(field, data):
    io = data.hierarchy.io
    if not data.NumberOfParticles > 0:
        return na.array([], dtype=dtype)
    try:
        return io._read_data_set(data, p_field).astype(dtype)
    except io._read_exception:
        pass
    # This is bad.  But it's the best idea I have right now.
    return data._read_data(p_field.replace("_"," ")).astype(dtype)

Convert Function Source

No source available.

particle_position_y

  • Particle Type: True

Field Source

def _Particles(field, data):
    io = data.hierarchy.io
    if not data.NumberOfParticles > 0:
        return na.array([], dtype=dtype)
    try:
        return io._read_data_set(data, p_field).astype(dtype)
    except io._read_exception:
        pass
    # This is bad.  But it's the best idea I have right now.
    return data._read_data(p_field.replace("_"," ")).astype(dtype)

Convert Function Source

No source available.

particle_position_z

  • Particle Type: True

Field Source

def _Particles(field, data):
    io = data.hierarchy.io
    if not data.NumberOfParticles > 0:
        return na.array([], dtype=dtype)
    try:
        return io._read_data_set(data, p_field).astype(dtype)
    except io._read_exception:
        pass
    # This is bad.  But it's the best idea I have right now.
    return data._read_data(p_field.replace("_"," ")).astype(dtype)

Convert Function Source

No source available.

particle_type

  • Particle Type: True

Field Source

def _Particles(field, data):
    io = data.hierarchy.io
    if not data.NumberOfParticles > 0:
        return na.array([], dtype=dtype)
    try:
        return io._read_data_set(data, p_field).astype(dtype)
    except io._read_exception:
        pass
    # This is bad.  But it's the best idea I have right now.
    return data._read_data(p_field.replace("_"," ")).astype(dtype)

Convert Function Source

No source available.

particle_velocity_x

  • Particle Type: True

Field Source

def _Particles(field, data):
    io = data.hierarchy.io
    if not data.NumberOfParticles > 0:
        return na.array([], dtype=dtype)
    try:
        return io._read_data_set(data, p_field).astype(dtype)
    except io._read_exception:
        pass
    # This is bad.  But it's the best idea I have right now.
    return data._read_data(p_field.replace("_"," ")).astype(dtype)

Convert Function Source

def _convert_p_vel(data):
    return data.convert("%s-velocity" % ax)

particle_velocity_y

  • Particle Type: True

Field Source

def _Particles(field, data):
    io = data.hierarchy.io
    if not data.NumberOfParticles > 0:
        return na.array([], dtype=dtype)
    try:
        return io._read_data_set(data, p_field).astype(dtype)
    except io._read_exception:
        pass
    # This is bad.  But it's the best idea I have right now.
    return data._read_data(p_field.replace("_"," ")).astype(dtype)

Convert Function Source

def _convert_p_vel(data):
    return data.convert("%s-velocity" % ax)

particle_velocity_z

  • Particle Type: True

Field Source

def _Particles(field, data):
    io = data.hierarchy.io
    if not data.NumberOfParticles > 0:
        return na.array([], dtype=dtype)
    try:
        return io._read_data_set(data, p_field).astype(dtype)
    except io._read_exception:
        pass
    # This is bad.  But it's the best idea I have right now.
    return data._read_data(p_field.replace("_"," ")).astype(dtype)

Convert Function Source

def _convert_p_vel(data):
    return data.convert("%s-velocity" % ax)

tempContours

  • Particle Type: False

Field Source

def _Contours(field, data):
    return na.ones(data["Density"].shape)*-1

Convert Function Source

No source available.

x

  • Particle Type: False

Field Source

def _coordX(field, data):
    dim = data.ActiveDimensions[0]
    return (na.ones(data.ActiveDimensions, dtype='float64')
                   * na.arange(data.ActiveDimensions[0])[:,None,None]
            +0.5) * data['dx'] + data.LeftEdge[0]

Convert Function Source

No source available.

y

  • Particle Type: False

Field Source

def _coordY(field, data):
    dim = data.ActiveDimensions[1]
    return (na.ones(data.ActiveDimensions, dtype='float64')
                   * na.arange(data.ActiveDimensions[1])[None,:,None]
            +0.5) * data['dy'] + data.LeftEdge[1]

Convert Function Source

No source available.

z

  • Particle Type: False

Field Source

def _coordZ(field, data):
    dim = data.ActiveDimensions[2]
    return (na.ones(data.ActiveDimensions, dtype='float64')
                   * na.arange(data.ActiveDimensions[2])[None,None,:]
            +0.5) * data['dz'] + data.LeftEdge[2]

Convert Function Source

No source available.

Enzo Field List

ComovingDensity

  • Units: \rm{g}/\rm{cm}^3
  • Particle Type: False

Field Source

def _ComovingDensity(field,data):
    ef = (1.0 + data.pf["CosmologyCurrentRedshift"])**3.0
    return data["Density"]/ef

Convert Function Source

No source available.

Comoving_DII_Density

  • Particle Type: False

Field Source

def _SpeciesComovingDensity(field, data):
    sp = field.name.split("_")[0] + "_Density"
    ef = (1.0 + data.pf["CosmologyCurrentRedshift"])**3.0
    return data[sp]/ef

Convert Function Source

No source available.

Comoving_DI_Density

  • Particle Type: False

Field Source

def _SpeciesComovingDensity(field, data):
    sp = field.name.split("_")[0] + "_Density"
    ef = (1.0 + data.pf["CosmologyCurrentRedshift"])**3.0
    return data[sp]/ef

Convert Function Source

No source available.

Comoving_Electron_Density

  • Particle Type: False

Field Source

def _SpeciesComovingDensity(field, data):
    sp = field.name.split("_")[0] + "_Density"
    ef = (1.0 + data.pf["CosmologyCurrentRedshift"])**3.0
    return data[sp]/ef

Convert Function Source

No source available.

Comoving_H2II_Density

  • Particle Type: False

Field Source

def _SpeciesComovingDensity(field, data):
    sp = field.name.split("_")[0] + "_Density"
    ef = (1.0 + data.pf["CosmologyCurrentRedshift"])**3.0
    return data[sp]/ef

Convert Function Source

No source available.

Comoving_H2I_Density

  • Particle Type: False

Field Source

def _SpeciesComovingDensity(field, data):
    sp = field.name.split("_")[0] + "_Density"
    ef = (1.0 + data.pf["CosmologyCurrentRedshift"])**3.0
    return data[sp]/ef

Convert Function Source

No source available.

Comoving_HDI_Density

  • Particle Type: False

Field Source

def _SpeciesComovingDensity(field, data):
    sp = field.name.split("_")[0] + "_Density"
    ef = (1.0 + data.pf["CosmologyCurrentRedshift"])**3.0
    return data[sp]/ef

Convert Function Source

No source available.

Comoving_HII_Density

  • Particle Type: False

Field Source

def _SpeciesComovingDensity(field, data):
    sp = field.name.split("_")[0] + "_Density"
    ef = (1.0 + data.pf["CosmologyCurrentRedshift"])**3.0
    return data[sp]/ef

Convert Function Source

No source available.

Comoving_HI_Density

  • Particle Type: False

Field Source

def _SpeciesComovingDensity(field, data):
    sp = field.name.split("_")[0] + "_Density"
    ef = (1.0 + data.pf["CosmologyCurrentRedshift"])**3.0
    return data[sp]/ef

Convert Function Source

No source available.

Comoving_HM_Density

  • Particle Type: False

Field Source

def _SpeciesComovingDensity(field, data):
    sp = field.name.split("_")[0] + "_Density"
    ef = (1.0 + data.pf["CosmologyCurrentRedshift"])**3.0
    return data[sp]/ef

Convert Function Source

No source available.

Comoving_HeIII_Density

  • Particle Type: False

Field Source

def _SpeciesComovingDensity(field, data):
    sp = field.name.split("_")[0] + "_Density"
    ef = (1.0 + data.pf["CosmologyCurrentRedshift"])**3.0
    return data[sp]/ef

Convert Function Source

No source available.

Comoving_HeII_Density

  • Particle Type: False

Field Source

def _SpeciesComovingDensity(field, data):
    sp = field.name.split("_")[0] + "_Density"
    ef = (1.0 + data.pf["CosmologyCurrentRedshift"])**3.0
    return data[sp]/ef

Convert Function Source

No source available.

Comoving_HeI_Density

  • Particle Type: False

Field Source

def _SpeciesComovingDensity(field, data):
    sp = field.name.split("_")[0] + "_Density"
    ef = (1.0 + data.pf["CosmologyCurrentRedshift"])**3.0
    return data[sp]/ef

Convert Function Source

No source available.

Comoving_Metal_Density

  • Particle Type: False

Field Source

def _SpeciesComovingDensity(field, data):
    sp = field.name.split("_")[0] + "_Density"
    ef = (1.0 + data.pf["CosmologyCurrentRedshift"])**3.0
    return data[sp]/ef

Convert Function Source

No source available.

Comoving_PreShock_Density

  • Particle Type: False

Field Source

def _SpeciesComovingDensity(field, data):
    sp = field.name.split("_")[0] + "_Density"
    ef = (1.0 + data.pf["CosmologyCurrentRedshift"])**3.0
    return data[sp]/ef

Convert Function Source

No source available.

Cooling_Time

  • Units: \rm{s}
  • Particle Type: False

Field Source

def _Cooling_Time(field, data):
    return data["Cooling_Time"]

Convert Function Source

No source available.

DII_Density

  • Units: \rm{g}/\rm{cm}^3
  • Projected Units: \rm{g}/\rm{cm}^2
  • Particle Type: False

Field Source

No source available.

DII_Fraction

  • Particle Type: False

Field Source

def _SpeciesFraction(field, data):
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/data["Density"]

Convert Function Source

No source available.

DII_NumberDensity

  • Particle Type: False

Field Source

def _SpeciesNumberDensity(field, data):
    species = field.name.split("_")[0]
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/_speciesMass[species]

Convert Function Source

def _ConvertNumberDensity(data):
    return 1.0/mh

DI_Density

  • Units: \rm{g}/\rm{cm}^3
  • Projected Units: \rm{g}/\rm{cm}^2
  • Particle Type: False

Field Source

No source available.

DI_Fraction

  • Particle Type: False

Field Source

def _SpeciesFraction(field, data):
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/data["Density"]

Convert Function Source

No source available.

DI_NumberDensity

  • Particle Type: False

Field Source

def _SpeciesNumberDensity(field, data):
    species = field.name.split("_")[0]
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/_speciesMass[species]

Convert Function Source

def _ConvertNumberDensity(data):
    return 1.0/mh

Dark_Matter_Density

  • Particle Type: False

Field Source

No source available.

Density

  • Units: \rm{g}/\rm{cm}^3
  • Projected Units: \rm{g}/\rm{cm}^2
  • Particle Type: False

Field Source

No source available.

Electron_Density

  • Units: \rm{g}/\rm{cm}^3
  • Projected Units: \rm{g}/\rm{cm}^2
  • Particle Type: False

Field Source

No source available.

Electron_Fraction

  • Particle Type: False

Field Source

def _SpeciesFraction(field, data):
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/data["Density"]

Convert Function Source

No source available.

Electron_NumberDensity

  • Particle Type: False

Field Source

def _SpeciesNumberDensity(field, data):
    species = field.name.split("_")[0]
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/_speciesMass[species]

Convert Function Source

def _ConvertNumberDensity(data):
    return 1.0/mh

GasEnergy

  • Units: \rm{ergs}/\rm{g}
  • Particle Type: False

Field Source

def _GasEnergy(field, data):
    return data["Gas_Energy"] / _convertEnergy(data)

Convert Function Source

def _convertEnergy(data):
    return data.convert("x-velocity")**2.0

Gas_Energy

  • Units: \rm{ergs}/\rm{g}
  • Particle Type: False

Field Source

def _Gas_Energy(field, data):
    return data["GasEnergy"] / _convertEnergy(data)

Convert Function Source

def _convertEnergy(data):
    return data.convert("x-velocity")**2.0

H2II_Density

  • Units: \rm{g}/\rm{cm}^3
  • Projected Units: \rm{g}/\rm{cm}^2
  • Particle Type: False

Field Source

No source available.

H2II_Fraction

  • Particle Type: False

Field Source

def _SpeciesFraction(field, data):
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/data["Density"]

Convert Function Source

No source available.

H2II_NumberDensity

  • Particle Type: False

Field Source

def _SpeciesNumberDensity(field, data):
    species = field.name.split("_")[0]
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/_speciesMass[species]

Convert Function Source

def _ConvertNumberDensity(data):
    return 1.0/mh

H2I_Density

  • Units: \rm{g}/\rm{cm}^3
  • Projected Units: \rm{g}/\rm{cm}^2
  • Particle Type: False

Field Source

No source available.

H2I_Fraction

  • Particle Type: False

Field Source

def _SpeciesFraction(field, data):
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/data["Density"]

Convert Function Source

No source available.

H2I_NumberDensity

  • Particle Type: False

Field Source

def _SpeciesNumberDensity(field, data):
    species = field.name.split("_")[0]
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/_speciesMass[species]

Convert Function Source

def _ConvertNumberDensity(data):
    return 1.0/mh

HDI_Density

  • Units: \rm{g}/\rm{cm}^3
  • Projected Units: \rm{g}/\rm{cm}^2
  • Particle Type: False

Field Source

No source available.

HDI_Fraction

  • Particle Type: False

Field Source

def _SpeciesFraction(field, data):
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/data["Density"]

Convert Function Source

No source available.

HDI_NumberDensity

  • Particle Type: False

Field Source

def _SpeciesNumberDensity(field, data):
    species = field.name.split("_")[0]
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/_speciesMass[species]

Convert Function Source

def _ConvertNumberDensity(data):
    return 1.0/mh

HII_Density

  • Units: \rm{g}/\rm{cm}^3
  • Projected Units: \rm{g}/\rm{cm}^2
  • Particle Type: False

Field Source

No source available.

HII_Fraction

  • Particle Type: False

Field Source

def _SpeciesFraction(field, data):
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/data["Density"]

Convert Function Source

No source available.

HII_NumberDensity

  • Particle Type: False

Field Source

def _SpeciesNumberDensity(field, data):
    species = field.name.split("_")[0]
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/_speciesMass[species]

Convert Function Source

def _ConvertNumberDensity(data):
    return 1.0/mh

HI_Density

  • Units: \rm{g}/\rm{cm}^3
  • Projected Units: \rm{g}/\rm{cm}^2
  • Particle Type: False

Field Source

No source available.

HI_Fraction

  • Particle Type: False

Field Source

def _SpeciesFraction(field, data):
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/data["Density"]

Convert Function Source

No source available.

HI_NumberDensity

  • Particle Type: False

Field Source

def _SpeciesNumberDensity(field, data):
    species = field.name.split("_")[0]
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/_speciesMass[species]

Convert Function Source

def _ConvertNumberDensity(data):
    return 1.0/mh

HM_Density

  • Units: \rm{g}/\rm{cm}^3
  • Projected Units: \rm{g}/\rm{cm}^2
  • Particle Type: False

Field Source

No source available.

HM_Fraction

  • Particle Type: False

Field Source

def _SpeciesFraction(field, data):
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/data["Density"]

Convert Function Source

No source available.

HM_NumberDensity

  • Particle Type: False

Field Source

def _SpeciesNumberDensity(field, data):
    species = field.name.split("_")[0]
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/_speciesMass[species]

Convert Function Source

def _ConvertNumberDensity(data):
    return 1.0/mh

HeIII_Density

  • Units: \rm{g}/\rm{cm}^3
  • Projected Units: \rm{g}/\rm{cm}^2
  • Particle Type: False

Field Source

No source available.

HeIII_Fraction

  • Particle Type: False

Field Source

def _SpeciesFraction(field, data):
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/data["Density"]

Convert Function Source

No source available.

HeIII_NumberDensity

  • Particle Type: False

Field Source

def _SpeciesNumberDensity(field, data):
    species = field.name.split("_")[0]
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/_speciesMass[species]

Convert Function Source

def _ConvertNumberDensity(data):
    return 1.0/mh

HeII_Density

  • Units: \rm{g}/\rm{cm}^3
  • Projected Units: \rm{g}/\rm{cm}^2
  • Particle Type: False

Field Source

No source available.

HeII_Fraction

  • Particle Type: False

Field Source

def _SpeciesFraction(field, data):
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/data["Density"]

Convert Function Source

No source available.

HeII_NumberDensity

  • Particle Type: False

Field Source

def _SpeciesNumberDensity(field, data):
    species = field.name.split("_")[0]
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/_speciesMass[species]

Convert Function Source

def _ConvertNumberDensity(data):
    return 1.0/mh

HeI_Density

  • Units: \rm{g}/\rm{cm}^3
  • Projected Units: \rm{g}/\rm{cm}^2
  • Particle Type: False

Field Source

No source available.

HeI_Fraction

  • Particle Type: False

Field Source

def _SpeciesFraction(field, data):
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/data["Density"]

Convert Function Source

No source available.

HeI_NumberDensity

  • Particle Type: False

Field Source

def _SpeciesNumberDensity(field, data):
    species = field.name.split("_")[0]
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/_speciesMass[species]

Convert Function Source

def _ConvertNumberDensity(data):
    return 1.0/mh

IsStarParticle

  • Particle Type: True

Field Source

def _IsStarParticle(field, data):
    is_star = (data['creation_time'] > 0).astype('float64')
    return is_star

Convert Function Source

No source available.

Metal_Density

  • Units: \rm{g}/\rm{cm}^3
  • Projected Units: \rm{g}/\rm{cm}^2
  • Particle Type: False

Field Source

No source available.

Metal_Fraction

  • Particle Type: False

Field Source

def _SpeciesFraction(field, data):
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/data["Density"]

Convert Function Source

No source available.

Metallicity

  • Units: Z_{\rm{\odot}}
  • Particle Type: False

Field Source

def _Metallicity(field, data):
    return data["Metal_Fraction"]

Convert Function Source

def _ConvertMetallicity(data):
    return 49.0196 # 1 / 0.0204

Metallicity3

  • Units: Z_{\rm{\odot}}
  • Particle Type: False

Field Source

def _Metallicity3(field, data):
    return data["SN_Colour"]/data["Density"]

Convert Function Source

def _ConvertMetallicity(data):
    return 49.0196 # 1 / 0.0204

NumberDensity

  • Units: \rm{cm}^{-3}
  • Particle Type: False

Field Source

def _NumberDensity(field, data):
    # We can assume that we at least have Density
    # We should actually be guaranteeing the presence of a .shape attribute,
    # but I am not currently implementing that
    fieldData = na.zeros(data["Density"].shape,
                         dtype = data["Density"].dtype)
    if data.pf["MultiSpecies"] == 0:
        if data.has_field_parameter("mu"):
            mu = data.get_field_parameter("mu")
        else:
            mu = 0.6
        fieldData += data["Density"] / mu
    if data.pf["MultiSpecies"] > 0:
        fieldData += data["HI_Density"] / 1.0
        fieldData += data["HII_Density"] / 1.0
        fieldData += data["HeI_Density"] / 4.0
        fieldData += data["HeII_Density"] / 4.0
        fieldData += data["HeIII_Density"] / 4.0
        fieldData += data["Electron_Density"] / 1.0
    if data.pf["MultiSpecies"] > 1:
        fieldData += data["HM_Density"] / 1.0
        fieldData += data["H2I_Density"] / 2.0
        fieldData += data["H2II_Density"] / 2.0
    if data.pf["MultiSpecies"] > 2:
        fieldData += data["DI_Density"] / 2.0
        fieldData += data["DII_Density"] / 2.0
        fieldData += data["HDI_Density"] / 3.0
    return fieldData

Convert Function Source

def _ConvertNumberDensity(data):
    return 1.0/mh

Overdensity

  • Particle Type: False

Field Source

def Overdensity(field,data):
    return (data['Density'] + data['Dark_Matter_Density']) / \
        (rho_crit_now * (data.pf['CosmologyHubbleConstantNow']**2) * ((1+data.pf['CosmologyCurrentRedshift'])**3))

Convert Function Source

No source available.

PreShock_Density

  • Units: \rm{g}/\rm{cm}^3
  • Projected Units: \rm{g}/\rm{cm}^2
  • Particle Type: False

Field Source

No source available.

PreShock_Fraction

  • Particle Type: False

Field Source

def _SpeciesFraction(field, data):
    sp = field.name.split("_")[0] + "_Density"
    return data[sp]/data["Density"]

Convert Function Source

No source available.

StarAgeYears

  • Units: \mathrm{yr}
  • Particle Type: False

Field Source

def _StarAge(field, data):
    star_age = na.zeros(data['StarCreationTimeYears'].shape)
    with_stars = data['StarCreationTimeYears'] > 0
    star_age[with_stars] = data.pf.time_units['years'] * \
        data.pf["InitialTime"] - \
        data['StarCreationTimeYears'][with_stars]
    return star_age

Convert Function Source

No source available.

StarCreationTimeYears

  • Units: \mathrm{yr}
  • Particle Type: False

Field Source

def _StarCreationTime(field, data):
    return data['star_creation_time']

Convert Function Source

def _ConvertEnzoTimeYears(data):
    return data.pf.time_units['years']

StarDynamicalTimeYears

  • Units: \mathrm{yr}
  • Particle Type: False

Field Source

def _StarDynamicalTime(field, data):
    return data['star_dynamical_time']

Convert Function Source

def _ConvertEnzoTimeYears(data):
    return data.pf.time_units['years']

StarMetallicity

  • Units: Z_{\rm{\odot}}
  • Particle Type: False

Field Source

def _StarMetallicity(field, data):
    return data['star_metallicity_fraction']

Convert Function Source

def _ConvertMetallicity(data):
    return 49.0196 # 1 / 0.0204

Temperature

  • Units: \rm{K}
  • Particle Type: False

Field Source

No source available.

ThermalEnergy

  • Units: \rm{ergs}/\rm{cm^3}
  • Particle Type: False

Field Source

def _ThermalEnergy(field, data):
    if data.pf["HydroMethod"] == 2:
        return data["Total_Energy"]
    else:
        if data.pf["DualEnergyFormalism"]:
            return data["GasEnergy"]
        else:
            return data["Total_Energy"] - 0.5*(
                   data["x-velocity"]**2.0
                 + data["y-velocity"]**2.0
                 + data["z-velocity"]**2.0 )

Convert Function Source

No source available.

TotalEnergy

  • Units: \rm{ergs}/\rm{g}
  • Particle Type: False

Field Source

def _TotalEnergy(field, data):
    return data["Total_Energy"] / _convertEnergy(data)

Convert Function Source

def _convertEnergy(data):
    return data.convert("x-velocity")**2.0

Total_Energy

  • Units: \rm{ergs}/\rm{g}
  • Particle Type: False

Field Source

def _Total_Energy(field, data):
    return data["TotalEnergy"] / _convertEnergy(data)

Convert Function Source

def _convertEnergy(data):
    return data.convert("x-velocity")**2.0

dm_density_pyx

  • Particle Type: False

Field Source

def _dmpdensity_pyx(field, data):
    blank = na.zeros(data.ActiveDimensions, dtype='float32')
    if data.NumberOfParticles == 0: return blank
    filter = data['creation_time'] <= 0.0
    if not filter.any(): return blank
    CICDeposit_3(data["particle_position_x"][filter].astype(na.float64),
                 data["particle_position_y"][filter].astype(na.float64),
                 data["particle_position_z"][filter].astype(na.float64),
                 data["particle_mass"][filter].astype(na.float32),
                 na.int64(na.where(filter)[0].size),
                 blank, na.array(data.LeftEdge).astype(na.float64),
                 na.array(data.ActiveDimensions).astype(na.int32),
                 na.float64(data['dx']))
    return blank

Convert Function Source

def _convertDensity(data):
    return data.convert("Density")

particle_density

  • Particle Type: False

Field Source

def _pdensity(field, data):
    blank = na.zeros(data.ActiveDimensions, dtype='float32', order="FORTRAN")
    if data.NumberOfParticles == 0: return blank
    cic_deposit.cic_deposit(data["particle_position_x"],
                            data["particle_position_y"],
                            data["particle_position_z"], 3,
                            data["particle_mass"],
                            blank, data.LeftEdge, data['dx'])
    return blank

Convert Function Source

def _convertDensity(data):
    return data.convert("Density")

star_creation_time

  • Particle Type: False

Field Source

def _star_field(field, data):
    """
    Create a grid field for star quantities, weighted by star mass.
    """
    particle_field = field.name[5:]
    top = na.zeros(data.ActiveDimensions, dtype='float32')
    if data.NumberOfParticles == 0: return top
    filter = data['creation_time'] > 0.0
    if not filter.any(): return top
    particle_field_data = data[particle_field][filter] * data['particle_mass'][filter]
    CICDeposit_3(data["particle_position_x"][filter].astype(na.float64),
                 data["particle_position_y"][filter].astype(na.float64),
                 data["particle_position_z"][filter].astype(na.float64),
                 particle_field_data.astype(na.float32),
                 na.int64(na.where(filter)[0].size),
                 top, na.array(data.LeftEdge).astype(na.float64),
                 na.array(data.ActiveDimensions).astype(na.int32),
                 na.float64(data['dx']))
    del particle_field_data

    bottom = na.zeros(data.ActiveDimensions, dtype='float32')
    CICDeposit_3(data["particle_position_x"][filter].astype(na.float64),
                 data["particle_position_y"][filter].astype(na.float64),
                 data["particle_position_z"][filter].astype(na.float64),
                 data["particle_mass"][filter].astype(na.float32),
                 na.int64(na.where(filter)[0].size),
                 bottom, na.array(data.LeftEdge).astype(na.float64),
                 na.array(data.ActiveDimensions).astype(na.int32),
                 na.float64(data['dx']))

    top[bottom == 0] = 0.0
    bnz = bottom.nonzero()
    top[bnz] /= bottom[bnz]
    return top

Convert Function Source

No source available.

star_density_pyx

  • Particle Type: False

Field Source

def _spdensity_pyx(field, data):
    blank = na.zeros(data.ActiveDimensions, dtype='float32')
    if data.NumberOfParticles == 0: return blank
    filter = data['creation_time'] > 0.0
    if not filter.any(): return blank
    CICDeposit_3(data["particle_position_x"][filter].astype(na.float64),
                 data["particle_position_y"][filter].astype(na.float64),
                 data["particle_position_z"][filter].astype(na.float64),
                 data["particle_mass"][filter].astype(na.float32),
                 na.int64(na.where(filter)[0].size),
                 blank, na.array(data.LeftEdge).astype(na.float64),
                 na.array(data.ActiveDimensions).astype(na.int32),
                 na.float64(data['dx']))
    return blank

Convert Function Source

def _convertDensity(data):
    return data.convert("Density")

star_dynamical_time

  • Particle Type: False

Field Source

def _star_field(field, data):
    """
    Create a grid field for star quantities, weighted by star mass.
    """
    particle_field = field.name[5:]
    top = na.zeros(data.ActiveDimensions, dtype='float32')
    if data.NumberOfParticles == 0: return top
    filter = data['creation_time'] > 0.0
    if not filter.any(): return top
    particle_field_data = data[particle_field][filter] * data['particle_mass'][filter]
    CICDeposit_3(data["particle_position_x"][filter].astype(na.float64),
                 data["particle_position_y"][filter].astype(na.float64),
                 data["particle_position_z"][filter].astype(na.float64),
                 particle_field_data.astype(na.float32),
                 na.int64(na.where(filter)[0].size),
                 top, na.array(data.LeftEdge).astype(na.float64),
                 na.array(data.ActiveDimensions).astype(na.int32),
                 na.float64(data['dx']))
    del particle_field_data

    bottom = na.zeros(data.ActiveDimensions, dtype='float32')
    CICDeposit_3(data["particle_position_x"][filter].astype(na.float64),
                 data["particle_position_y"][filter].astype(na.float64),
                 data["particle_position_z"][filter].astype(na.float64),
                 data["particle_mass"][filter].astype(na.float32),
                 na.int64(na.where(filter)[0].size),
                 bottom, na.array(data.LeftEdge).astype(na.float64),
                 na.array(data.ActiveDimensions).astype(na.int32),
                 na.float64(data['dx']))

    top[bottom == 0] = 0.0
    bnz = bottom.nonzero()
    top[bnz] /= bottom[bnz]
    return top

Convert Function Source

No source available.

star_metallicity_fraction

  • Particle Type: False

Field Source

def _star_field(field, data):
    """
    Create a grid field for star quantities, weighted by star mass.
    """
    particle_field = field.name[5:]
    top = na.zeros(data.ActiveDimensions, dtype='float32')
    if data.NumberOfParticles == 0: return top
    filter = data['creation_time'] > 0.0
    if not filter.any(): return top
    particle_field_data = data[particle_field][filter] * data['particle_mass'][filter]
    CICDeposit_3(data["particle_position_x"][filter].astype(na.float64),
                 data["particle_position_y"][filter].astype(na.float64),
                 data["particle_position_z"][filter].astype(na.float64),
                 particle_field_data.astype(na.float32),
                 na.int64(na.where(filter)[0].size),
                 top, na.array(data.LeftEdge).astype(na.float64),
                 na.array(data.ActiveDimensions).astype(na.int32),
                 na.float64(data['dx']))
    del particle_field_data

    bottom = na.zeros(data.ActiveDimensions, dtype='float32')
    CICDeposit_3(data["particle_position_x"][filter].astype(na.float64),
                 data["particle_position_y"][filter].astype(na.float64),
                 data["particle_position_z"][filter].astype(na.float64),
                 data["particle_mass"][filter].astype(na.float32),
                 na.int64(na.where(filter)[0].size),
                 bottom, na.array(data.LeftEdge).astype(na.float64),
                 na.array(data.ActiveDimensions).astype(na.int32),
                 na.float64(data['dx']))

    top[bottom == 0] = 0.0
    bnz = bottom.nonzero()
    top[bnz] /= bottom[bnz]
    return top

Convert Function Source

No source available.

x-momentum

  • Units: \rm{g}/\rm{cm}^3
  • Particle Type: False

Field Source

No source available.

x-velocity

  • Units: \rm{cm}/\rm{s}
  • Particle Type: False

Field Source

No source available.

y-momentum

  • Units: \rm{g}/\rm{cm}^3
  • Particle Type: False

Field Source

No source available.

y-velocity

  • Units: \rm{cm}/\rm{s}
  • Particle Type: False

Field Source

No source available.

z-momentum

  • Units: \rm{g}/\rm{cm}^3
  • Particle Type: False

Field Source

No source available.

z-velocity

  • Units: \rm{cm}/\rm{s}
  • Particle Type: False

Field Source

No source available.

Orion Field List

Density

  • Particle Type: False

Field Source

No source available.

Pressure

  • Units: \rm{dyne}/\rm{cm}^{2}
  • Particle Type: False

Field Source

def _Pressure(field,data):
    """M{(Gamma-1.0)*e, where e is thermal energy density
       NB: this will need to be modified for radiation
    """
    return (data.pf["Gamma"] - 1.0)*data["ThermalEnergy"]

Convert Function Source

No source available.

Temperature

  • Units: \rm{Kelvin}
  • Particle Type: False

Field Source

def _Temperature(field,data):
    return (data.pf["Gamma"]-1.0)*data.pf["mu"]*mh*data["ThermalEnergy"]/(kboltz*data["Density"])

Convert Function Source

No source available.

ThermalEnergy

  • Units: \rm{ergs}/\rm{cm^3}
  • Particle Type: False

Field Source

def _ThermalEnergy(field, data):
    """generate thermal (gas energy). Dual Energy Formalism was
        implemented by Stella, but this isn't how it's called, so I'll
        leave that commented out for now.
    """
    #if data.pf["DualEnergyFormalism"]:
    #    return data["Gas_Energy"]
    #else:
    return data["Total_Energy"] - 0.5 * data["density"] * (
        data["x-velocity"]**2.0
        + data["y-velocity"]**2.0
        + data["z-velocity"]**2.0 )

Convert Function Source

No source available.

Total_Energy

  • Particle Type: False

Field Source

No source available.

density

  • Units: \rm{g}/\rm{cm}^3
  • Projected Units: \rm{g}/\rm{cm}^2
  • Particle Type: False

Field Source

No source available.

eden

  • Units: \rm{erg}/\rm{cm}^3
  • Particle Type: False

Field Source

No source available.

temperature

  • Particle Type: False

Field Source

No source available.

x-momentum

  • Particle Type: False

Field Source

No source available.

x-velocity

  • Units: \rm{cm}/\rm{s}
  • Particle Type: False

Field Source

def _xVelocity(field, data):
    """generate x-velocity from x-momentum and density

    """
    return data["xmom"]/data["density"]

Convert Function Source

No source available.

xmom

  • Units: \rm{g}/\rm{cm^2\ s}
  • Particle Type: False

Field Source

No source available.

xvel

  • Particle Type: False

Field Source

No source available.

y-momentum

  • Particle Type: False

Field Source

No source available.

y-velocity

  • Units: \rm{cm}/\rm{s}
  • Particle Type: False

Field Source

def _yVelocity(field,data):
    """generate y-velocity from y-momentum and density

    """
    #try:
    #    return data["xvel"]
    #except KeyError:
    return data["ymom"]/data["density"]

Convert Function Source

No source available.

ymom

  • Units: \rm{gm}/\rm{cm^2\ s}
  • Particle Type: False

Field Source

No source available.

yvel

  • Particle Type: False

Field Source

No source available.

z-momentum

  • Particle Type: False

Field Source

No source available.

z-velocity

  • Units: \rm{cm}/\rm{s}
  • Particle Type: False

Field Source

def _zVelocity(field,data):
    """generate z-velocity from z-momentum and density

    """
    return data["zmom"]/data["density"]

Convert Function Source

No source available.

zmom

  • Units: \rm{g}/\rm{cm^2\ s}
  • Particle Type: False

Field Source

No source available.

zvel

  • Particle Type: False

Field Source

No source available.

Table Of Contents

Previous topic

yt.raven.get_multi_plot

Next topic

FAQ

This Page