In yt 3.0, we want to make it easier to access "raw" simulation data that a code writes directly to disk. The new unit system makes it much simpler to convert back and forth between physical coordinates and the unscaled "raw" coordinate system used internally in the simulation code. In some cases, this conversion involves transforming to comoving coordinates, so that is also covered here.
Let's take a look at a cosmological enzo dataset to play with converting between physical units and code units:
import yt
ds = yt.load('Enzo_64/DD0043/data0043')
The conversion factors between Enzo's internal unit system and the physical CGS system are stored in the dataset's unit_registry
object. Code units have names like code_length
and code_time
. Let's take a look at the names of all of the code units, along with their CGS conversion factors for this cosmological enzo dataset:
reg = ds.unit_registry
for un in reg.keys():
if un.startswith('code_'):
fmt_tup = (un, reg.lut[un][0], str(reg.lut[un][1]))
print ("Unit name: {:<15}\nCGS conversion: {:<15}\nDimensions: {:<15}\n".format(*fmt_tup))
fmt_tup
Most of the time you will not have to deal with the unit registry. For example, the conversion factors to code units are stored as attributes of the dataset object:
print ("Length unit: ", ds.length_unit)
print ("Time unit: ", ds.time_unit)
print ("Mass unit: ", ds.mass_unit)
print ("Velocity unit: ", ds.velocity_unit)
Conversion factors will be supplied in CGS by default. We can also ask what the conversion factors are in code units.
print ("Length unit: ", ds.length_unit.in_units('code_length'))
print ("Time unit: ", ds.time_unit.in_units('code_time'))
print ("Mass unit: ", ds.mass_unit.in_units('code_mass'))
print ("Velocity unit: ", ds.velocity_unit.in_units('code_velocity'))
as expected, all the conversion factors are unity in code units.
We can also play with unit conversions on ds.domain_width
. First, we see for enzo how code length units are defined relative to the domain width:
ds.domain_width
ds.domain_width.in_cgs()
ds.domain_width.in_units('Mpccm/h')
This last example uses a cosmological unit. In english, I asked for the domain width in comoving megaparsecs, scaled as if the hubble constant were 100 km/s/Mpc. Although $h$ isn't really a unit, yt treats it as one for the purposes of the unit system.
As an aside, Darren Croton's research note on the history, use, and interpretation of $h$ as it appears in the astronomical literature is pretty much required reading for anyone who has to deal with factors of $h$ every now and then.
In yt, comoving length unit symbols are named following the pattern (length symbol)cm
, i.e. pccm
for comoving parsec or mcm
for a comoving meter. A comoving length unit is different from the normal length unit by a factor of $(1+z)$:
z = ds.current_redshift
print (ds.quan(1, 'Mpc')/ds.quan(1, 'Mpccm'))
print (1+z)
As we saw before, $h$ is treated like any other unit symbol. It has dimensionless
units, just like a scalar:
print (ds.quan(1, 'Mpc')/ds.quan(1, 'Mpc/h'))
print (ds.hubble_constant)
These units can be used in readily used in plots and anywhere a length unit is appropriate in yt.
slc = yt.SlicePlot(ds, 0, 'density', width=(128, 'Mpccm/h'))
slc.set_figure_size(6)
When you create a YTArray
without referring to a unit registry, yt uses the default unit registry, which does not include code units or comoving units.
from yt import YTQuantity
a = YTQuantity(3, 'cm')
print (a.units.registry.keys())
When a dataset is loaded, yt infers conversion factors from the internal simulation unit system to the CGS unit system. These conversion factors are stored in a unit_registry
along with conversion factors to the other known unit symbols. For the cosmological Enzo dataset we loaded earlier, we can see there are a number of additional unit symbols not defined in the default unit lookup table:
print (sorted([k for k in ds.unit_registry.keys() if k not in a.units.registry.keys()]))
Since code units do not appear in the default unit symbol lookup table, one must explicitly refer to a unit registry when creating a YTArray
to be able to convert to the unit system of a simulation.
To make this as clean as possible, there are array and quantity-creating convenience functions attached to the Dataset
object:
ds.arr()
ds.quan()
These functions make it straightforward to create arrays and quantities that can be converted to code units or comoving units. For example:
a = ds.quan(3, 'code_length')
print (a)
print (a.in_cgs())
print (a.in_units('Mpccm/h'))
b = ds.arr([3, 4, 5], 'Mpccm/h')
print (b)
print (b.in_cgs())
On occasion, you might have a dataset for a supported frontend that does not have the conversions to code units accessible (for example, Athena data) or you may want to change them outright. yt
provides a mechanism so that one may provide their own code unit definitions to load
, which override the default rules for a given frontend for defining code units. This is provided through the units_override
dictionary. We'll use an example of an Athena dataset. First, a call to load
without units_override
:
ds1 = yt.load("MHDSloshing/virgo_low_res.0054.vtk")
print (ds1.length_unit)
print (ds1.mass_unit)
print (ds1.time_unit)
sp1 = ds1.sphere("c",(0.1,"unitary"))
print (sp1["density"])
This is a galaxy cluster dataset, so it is not likely that the units of density are correct. We happen to know that the unit definitions are different, so we can override the units:
units_override = {"length_unit":(1.0,"Mpc"),
"time_unit":(1.0,"Myr"),
"mass_unit":(1.0e14,"Msun")}
units_override
can take the following keys:
length_unit
time_unit
mass_unit
magnetic_unit
temperature_unit
and the associated values can be (value, unit) tuples, YTQuantities
, or floats (in the latter case they are assumed to have the corresponding cgs unit).
ds2 = yt.load("MHDSloshing/virgo_low_res.0054.vtk", units_override=units_override)
print (ds2.length_unit)
print (ds2.mass_unit)
print (ds2.time_unit)
sp2 = ds2.sphere("c",(0.1,"unitary"))
print (sp2["density"])
This option should be used very carefully, and only if you know that the dataset does not provide units or that the unit definitions generated are incorrect for some reason.
(3)_Comoving_units_and_code_units.ipynb; 3)_Comoving_units_and_code_units_evaluated.ipynb; 3)_Comoving_units_and_code_units.py)
yt has additional capabilities to handle the comoving coordinate system used internally in cosmological simulations. Simulations that use comoving coordinates, all length units have three other counterparts corresponding to comoving units, scaled comoving units, and scaled proper units. In all cases ‘scaled’ units refer to scaling by the reduced Hubble parameter - i.e. the length unit is what it would be in a universe where Hubble’s parameter is 100 km/s/Mpc.
To access these different units, yt has a common naming system. Scaled units are denoted by
dividing by the scaled Hubble parameter h
(which is in itself a unit). Comoving
units are denoted by appending cm
to the end of the unit name.
Using the parsec as an example,
pc
pccm
pccm/h
pc/h