Symbolic units: yt.unitsΒΆ

Notebook

Dimensional analysis

The fastest way to get into the unit system is to explore the quantities that live in the yt.units namespace:

In [1]:
from yt.units import meter, gram, kilogram, second, joule
print (kilogram*meter**2/second**2 == joule)
print (kilogram*meter**2/second**2)
True
1.0 kg*m**2/s**2
In [2]:
from yt.units import m, kg, s, W
kg*m**2/s**3 == W
Out[2]:
array(True, dtype=bool)
In [3]:
from yt.units import kilometer
three_kilometers = 3*kilometer
print (three_kilometers)
3.0 km
In [4]:
from yt.units import gram, kilogram
print (gram+kilogram)

print (kilogram+gram)

print (kilogram/gram)
1001.0 g
1.001 kg
1000.0 dimensionless

These unit symbols are all instances of a new class we've added to yt 3.0, YTQuantity. YTQuantity is useful for storing a single data point.

In [5]:
type(kilogram)
Out[5]:
yt.units.yt_array.YTQuantity

We also provide YTArray, which can store arrays of quantities:

In [6]:
arr = [3,4,5]*kilogram

print (arr)

print (type(arr))
[ 3.  4.  5.] kg
<class 'yt.units.yt_array.YTArray'>

Creating arrays and quantities

Most people will interact with the new unit system using YTArray and YTQuantity. These are both subclasses of numpy's fast array type, ndarray, and can be used interchangably with other NumPy arrays. These new classes make use of the unit system to append unit metadata to the underlying ndarray. YTArray is intended to store array data, while YTQuantitity is intended to store scalars in a particular unit system.

There are two ways to create arrays and quantities. The first is to explicitly create it by calling the class constructor and supplying a unit string:

In [7]:
from yt.units.yt_array import YTArray

sample_array = YTArray([1,2,3], 'g/cm**3')

print (sample_array)
[ 1.  2.  3.] g/cm**3

The unit string can be an arbitrary combination of metric unit names. Just a few examples:

In [8]:
from yt.units.yt_array import YTQuantity
from yt.units import kboltz
from numpy.random import random
import numpy as np

print ("Length:")
print (YTQuantity(random(), 'm'))
print (YTQuantity(random(), 'cm'))
print (YTQuantity(random(), 'Mpc'))
print (YTQuantity(random(), 'AU'))
print ('')

print ("Time:")
print (YTQuantity(random(), 's'))
print (YTQuantity(random(), 'min'))
print (YTQuantity(random(), 'hr'))
print (YTQuantity(random(), 'day'))
print (YTQuantity(random(), 'yr'))
print ('')

print ("Mass:")
print (YTQuantity(random(), 'g'))
print (YTQuantity(random(), 'kg'))
print (YTQuantity(random(), 'Msun'))
print ('')

print ("Energy:")
print (YTQuantity(random(), 'erg'))
print (YTQuantity(random(), 'g*cm**2/s**2'))
print (YTQuantity(random(), 'eV'))
print (YTQuantity(random(), 'J'))
print ('')

print ("Temperature:")
print (YTQuantity(random(), 'K'))
print ((YTQuantity(random(), 'eV')/kboltz).in_cgs())
Length:
0.12985471124564796 m
0.22377174326354576 cm
0.8212243679210592 Mpc
0.13565541924712898 AU

Time:
0.25785985170790315 s
0.3403006842607613 min
0.9289241211122439 hr
0.5262333421539168 day
0.33882712911229396 yr

Mass:
0.10528700735352414 g
0.6483752053287599 kg
0.21827399006571246 Msun

Energy:
0.1604208827550596 erg
0.046548079143261756 cm**2*g/s**2
0.6559941387427474 eV
0.2389725922718804 J

Temperature:
0.9979997368329508 K
5919.783636928663 K

Dimensional arrays and quantities can also be created by multiplication with another array or quantity:

In [9]:
from yt.units import kilometer
print (kilometer)
1.0 km
In [10]:
three_kilometers = 3*kilometer
print (three_kilometers)
3.0 km

When working with a YTArray with complicated units, you can use unit_array and unit_quantity to conveniently apply units to data:

In [11]:
test_array = YTArray(np.random.random(20), 'erg/s')

print (test_array)
[ 0.77271739  0.7264211   0.25838029  0.90680308  0.49870713  0.35443695
  0.05498492  0.1130506   0.5386724   0.90003513  0.25478248  0.82998205
  0.84979257  0.83956161  0.15618761  0.70062988  0.13111743  0.32615018
  0.81905343  0.03367022] erg/s

unit_quantity returns a YTQuantity with a value of 1.0 and the same units as the array it is a attached to.

In [12]:
print (test_array.unit_quantity)
1.0 erg/s

unit_array returns a YTArray with the same units and shape as the array it is a attached to and with all values set to 1.0.

In [13]:
print (test_array.unit_array)
[ 1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.
  1.  1.] erg/s

These are useful when doing arithmetic:

In [14]:
print (test_array + 1.0*test_array.unit_quantity)
[ 1.77271739  1.7264211   1.25838029  1.90680308  1.49870713  1.35443695
  1.05498492  1.1130506   1.5386724   1.90003513  1.25478248  1.82998205
  1.84979257  1.83956161  1.15618761  1.70062988  1.13111743  1.32615018
  1.81905343  1.03367022] erg/s
In [15]:
print (test_array + np.arange(20)*test_array.unit_array)
[  0.77271739   1.7264211    2.25838029   3.90680308   4.49870713
   5.35443695   6.05498492   7.1130506    8.5386724    9.90003513
  10.25478248  11.82998205  12.84979257  13.83956161  14.15618761
  15.70062988  16.13111743  17.32615018  18.81905343  19.03367022] erg/s

For convenience, unit_quantity is also available via uq and unit_array is available via ua. You can use these arrays to create dummy arrays with the same units as another array - this is sometimes easier than manually creating a new array or quantity.

In [16]:
print (test_array.uq)

print (test_array.unit_quantity == test_array.uq)
1.0 erg/s
True
In [17]:
from numpy import array_equal

print (test_array.ua)

print (array_equal(test_array.ua, test_array.unit_array))
[ 1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.  1.
  1.  1.] erg/s
True

Unit metadata is encoded in the units attribute that hangs off of YTArray or YTQuantity instances:

In [18]:
from yt.units import kilometer, erg

print ("kilometer's units:", kilometer.units)
print ("kilometer's dimensions:", kilometer.units.dimensions)

print ('')

print ("erg's units:", erg.units)
print ("erg's dimensions: ", erg.units.dimensions)
kilometer's units: km
kilometer's dimensions: (length)

erg's units: erg
erg's dimensions:  (length)**2*(mass)/(time)**2

Arithmetic with YTQuantity and YTArray

Of course it wouldn't be very useful if all we could do is create data with units. The real power of the new unit system is that we can add, subtract, mutliply, and divide using quantities and dimensional arrays:

In [19]:
a = YTQuantity(3, 'cm')
b = YTQuantity(3, 'm')

print (a+b)
print (b+a)
print ('')

print ((a+b).in_units('ft'))
303.0 cm
3.03 m

9.940944881889763 ft
In [20]:
a = YTQuantity(42, 'mm')
b = YTQuantity(1, 's')

print (a/b)
print ((a/b).in_cgs())
print ((a/b).in_mks())
print ((a/b).in_units('km/s'))
print ('')

print (a*b)
print ((a*b).in_cgs())
print ((a*b).in_mks())
42.0 mm/s
4.2 cm/s
0.042 m/s
4.2e-05 km/s

42.0 mm*s
4.2 cm*s
0.042 m*s
In [21]:
m = YTQuantity(35, 'g')
a = YTQuantity(9.8, 'm/s**2')

print (m*a)
print ((m*a).in_units('dyne'))
343.0 g*m/s**2
34300.0 dyne
In [22]:
from yt.units import G, kboltz

print ("Newton's constant: ", G)
print ("Newton's constant in MKS: ", G.in_mks(), "\n")

print ("Boltzmann constant: ", kboltz)
print ("Boltzmann constant in MKS: ", kboltz.in_mks())
Newton's constant:  6.67384e-08 cm**3/(g*s**2)
Newton's constant in MKS:  6.67384e-11 m**3/(kg*s**2) 

Boltzmann constant:  1.3806488e-16 erg/K
Boltzmann constant in MKS:  1.3806488e-23 kg*m**2/(K*s**2)
In [23]:
rho = YTQuantity(1, 'g/cm**3')
t_ff = (G*rho)**(-0.5)

print (t_ff)
3870.901361178502 s

An exception is raised if we try to do a unit operation that doesn't make any sense:

In [24]:
from yt.utilities.exceptions import YTUnitOperationError

a = YTQuantity(3, 'm')
b = YTQuantity(5, 'erg')

try:
    print (a+b)
except YTUnitOperationError as e:
    print (e)
The addition operator for YTArrays with units (m) and (erg) is not well defined.

A plain ndarray or a YTArray created with empty units is treated as a dimensionless quantity and can be used in situations where unit consistency allows it to be used:

In [25]:
a = YTArray([1.,2.,3.], 'm')
b = np.array([2.,2.,2.])

print ("a:   ", a)
print ("b:   ", b)
print ("a*b: ", a*b)
a:    [ 1.  2.  3.] m
b:    [ 2.  2.  2.]
a*b:  [ 2.  4.  6.] m
In [26]:
c = YTArray([2,2,2])

print ("c:    ", c)
print ("a*c:  ", a*c)
c:     [ 2.  2.  2.] dimensionless
a*c:   [ 2.  4.  6.] m

Saving and Loading YTArrays to/from disk

YTArrays can be written to disk, to be loaded again to be used in yt or in a different context later. There are two formats that can be written to/read from: HDF5 and ASCII.

HDF5

To write to HDF5, use write_hdf5:

In [27]:
my_dens = YTArray(np.random.random(10), 'Msun/kpc**3')
my_temp = YTArray(np.random.random(10), 'K')
my_dens.write_hdf5("my_data.h5", dataset_name="density")
my_temp.write_hdf5("my_data.h5", dataset_name="temperature")

Where we used the dataset_name keyword argument to create a separate dataset for each array in the same file.

We can use the from_hdf5 classmethod to read the data back in:

In [28]:
read_dens = YTArray.from_hdf5("my_data.h5", dataset_name="density")
print (read_dens)
print (my_dens)
[ 0.97864639  0.11283437  0.9424725   0.1600274   0.12931168  0.09287471
  0.93297234  0.1214811   0.28221973  0.28219938] Msun/kpc**3
[ 0.97864639  0.11283437  0.9424725   0.1600274   0.12931168  0.09287471
  0.93297234  0.1214811   0.28221973  0.28219938] Msun/kpc**3

We can use the info keyword argument to write_hdf5 to write some additional data to the file, which will be stored as attributes of the dataset:

In [29]:
my_vels = YTArray(np.random.normal(10), 'km/s')
info = {"source":"galaxy cluster","user":"jzuhone"}
my_vels.write_hdf5("my_data.h5", dataset_name="velocity", info=info)

If you want to read/write a dataset from/to a specific group within the HDF5 file, use the group_name keyword:

In [30]:
my_vels.write_hdf5("data_in_group.h5", dataset_name="velocity", info=info, group_name="/data/fields")

where we have used the standard HDF5 slash notation for writing a group hierarchy (e.g., group within a group):

ASCII

To write one or more YTArrays to an ASCII text file, use yt.savetxt, which works a lot like NumPy's savetxt, except with units:

In [31]:
import yt
a = YTArray(np.random.random(size=10), "cm")
b = YTArray(np.random.random(size=10), "g")
c = YTArray(np.random.random(size=10), "s")
yt.savetxt("my_data.dat", [a,b,c], header='My cool data', footer='Data is over', delimiter="\t")

The file we wrote can then be easily used in other contexts, such as plotting in Gnuplot, or loading into a spreadsheet, or just for causal examination. We can quickly check it here:

In [32]:
%%bash 
more my_data.dat
::::::::::::::
my_data.dat
::::::::::::::
#My cool data
# Units
# cm	g	s
1.790070709828877060e-01	3.617349392898737692e-01	1.935540933902315519e-01
8.816076917726620721e-01	3.537521978119505528e-01	4.000387451314395548e-01
6.363512309271850409e-01	2.264057151788112510e-02	3.268207651487775589e-01
3.726790369651693524e-01	7.845609216083510029e-01	4.468130647047965365e-01
9.608941170864268022e-01	9.476803411703489388e-01	4.380286177888276233e-01
3.046693262522013335e-01	5.263235040333804626e-01	7.732731526107573528e-01
1.650222994083557770e-01	8.362779967355660204e-02	3.224219408887036842e-01
1.936902520762832092e-01	3.361878413910099894e-01	9.844138326884591672e-01
1.333276913387991813e-01	7.775933322579353657e-01	7.955059401022968046e-02
4.252518129006717107e-01	2.892294931626900878e-01	6.641053540807443367e-01
#Data is over

You can see that the header comes first, and then right before the data we have a subheader marking the units of each column. The footer comes after the data.

yt.loadtxt can be used to read the same data with units back in, or read data that has been generated from some other source. Just make sure it's in the format above. loadtxt can also selectively read from particular columns in the file with the usecols keyword argument:

In [33]:
bb, cc = yt.loadtxt("my_data.dat", usecols=(1,2), delimiter="\t")
print (bb)
print (b)
print ('')
print (cc)
print (c)
[ 0.36173494  0.3537522   0.02264057  0.78456092  0.94768034  0.5263235
  0.0836278   0.33618784  0.77759333  0.28922949] g
[ 0.36173494  0.3537522   0.02264057  0.78456092  0.94768034  0.5263235
  0.0836278   0.33618784  0.77759333  0.28922949] g

[ 0.19355409  0.40003875  0.32682077  0.44681306  0.43802862  0.77327315
  0.32242194  0.98441383  0.07955059  0.66410535] s
[ 0.19355409  0.40003875  0.32682077  0.44681306  0.43802862  0.77327315
  0.32242194  0.98441383  0.07955059  0.66410535] s

(1)_Symbolic_Units.ipynb; 1)_Symbolic_Units_evaluated.ipynb; 1)_Symbolic_Units.py)