Essential Framework Modules

array

Defines a shallow wrapper around numpy.ndarray for extra functionality like unit-tracking.

For most purposes, the differences between numpy.ndarray and array.SimArray are not important. However, when units are specified (by setting the units attribute), the behaviour is slightly different. In particular,

• it becomes impossible to add or subtract arrays with incompatible dimensions

>>> SimArray([1,2], "Mpc") + SimArray([1,2], "Msol"))
ValueError

• addition or subtraction causes auto-unit conversion. For example

>>> SimArray([1,2], "Mpc") + SimArray([1,2], "kpc")
SimArray([1.001, 1.002], "Mpc")

• Note that in this context the left value takes precedence in specifying the return units, so that reversing the order of the operation here would return results in kpc.

• If only one of the arrays specifies a Unit, no checking occurs and the unit of the returned array is assumed to be the same as the one specified input unit.

• Powers to single integer or rational powers will maintain unit tracking. Powers to float or other powers will not be able to do so.

>>> SimArray([1,2],"Msol Mpc**-3")**2
SimArray([1, 4], 'Msol**2 Mpc**-6')
>>> SimArray([1,2],"Msol Mpc**-3")**(1,3)
SimArray([ 1.,1.26], 'Msol**1/3 Mpc**-1')

Syntax above mirrors syntax in units module, where a length-two tuple can represent a rational number, in this case one third.

>>> SimArray([1.,2], "Msol Mpc**-3")**0.333
SimArray([ 1.,1.26])  # Lost track of units

Getting the array in specified units

Given an array, you can convert it in-place into units of your own chosing:

>>> x = SimArray([1,2], "Msol")
>>> x.convert_units('kg')
>>> print x
SimArray([  1.99e+30,   3.98e+30], 'kg')

Or you can leave the original array alone and get a copy in different units, correctly converted:

>>> x = SimArray([1,2], "Msol")
>>> print x.in_units("kg")
SimArray([  1.99e+30,   3.98e+30], 'kg')
>>> print x
SimArray([1,2], "Msol")

If the SimArray was created by a SimSnap (which is most likely), it has a pointer into the SimSnap’s properties so that the cosmological context is automatically fetched. For example, comoving -> physical conversions are correctly achieved:

>>> f = pynbody.load("fname")
>>> f['pos']
SimArray([[ 0.05419805, -0.0646539 , -0.15700017],
         [ 0.05169899, -0.06193341, -0.14475258],
         [ 0.05672406, -0.06384531, -0.15909944],
         ...,
         [ 0.0723075 , -0.07650762, -0.07657281],
         [ 0.07166634, -0.07453796, -0.08020873],
         [ 0.07165282, -0.07468577, -0.08020587]], '2.86e+04 kpc a')

>>> f['pos'].convert_units('kpc')
>>> f['pos']
SimArray([[ 1548.51403101, -1847.2525312 , -4485.71463308],
         [ 1477.1124212 , -1769.52421398, -4135.78377699],
         [ 1620.68592366, -1824.15000686, -4545.69387564],
         ...,
         [ 2065.9264273 , -2185.92982874, -2187.79225915],
         [ 2047.60759667, -2129.6537339 , -2291.6758134 ],
         [ 2047.2214441 , -2133.87693163, -2291.59406997]], 'kpc')

Specifying rules for ufunc’s

In general, it’s not possible to infer what the output units from a given ufunc should be. While numpy built-in ufuncs should be handled OK, other ufuncs will need their output units defined (otherwise a numpy.ndarray will be returned instead of our custom type.)

To do this, decorate a function with SimArray.ufunc_rule(ufunc). The function you define should take the same number of parameters as the ufunc. These will be the input parameters of the ufunc. You should return the correct units for the output, or raise units.UnitsException (in the latter case, the return array will be made into a numpy.ndarray.)

For example, here is the code for the correct addition/subtraction handler:

@SimArray.ufunc_rule(np.add)
@SimArray.ufunc_rule(np.subtract)
def _consistent_units(a,b) :

    # This will be called whenever the standard numpy ufuncs np.add
    # or np.subtract are called with parameters a,b.

    # You should always be ready for the inputs to have no units.

    a_units = a.units if hasattr(a, 'units') else None
    b_units = b.units if hasattr(b, 'units') else None

    # Now do the logic. If we're adding incompatible units,
    # we want just to get a plain numpy array out. If we only
    # know the units of one of the arrays, we assume the output
    # is in those units.

    if a_units is not None and b_units is not None :
        if a_units==b_units :
            return a_units
        else :
            raise units.UnitsException("Incompatible units")

    elif a_units is not None :
        return a_units
    else :
        return b_units

class pynbody.array.SimArray(data, units=None, sim=None, **kwargs)

Defines a shallow wrapper around numpy.ndarray for extra functionality like unit-tracking.

Attributes

T

The transposed array.

ancestor

Provides the basemost SimArray that an IndexedSimArray is based on.

base

Base object if memory is from some other object.

ctypes

An object to simplify the interaction of the array with the ctypes module.

data

Python buffer object pointing to the start of the array’s data.

derived

dtype

Data-type of the array’s elements.

family

flags

Information about the memory layout of the array.

flat

A 1-D iterator over the array.

imag

The imaginary part of the array.

itemsize

Length of one array element in bytes.

name

nbytes

Total bytes consumed by the elements of the array.

ndim

Number of array dimensions.

real

The real part of the array.

shape

Tuple of array dimensions.

sim

size

Number of elements in the array.

strides

Tuple of bytes to step in each dimension when traversing an array.

units

Methods

all([axis, out, keepdims])	Returns True if all elements evaluate to True.
any([axis, out, keepdims])	Returns True if any of the elements of a evaluate to True.
argmax([axis, out])	Return indices of the maximum values along the given axis.
argmin([axis, out])	Return indices of the minimum values along the given axis of a.
argpartition(kth[, axis, kind, order])	Returns the indices that would partition this array.
argsort([axis, kind, order])	Returns the indices that would sort this array.
astype(dtype[, order, casting, subok, copy])	Copy of the array, cast to a specified type.
byteswap([inplace])	Swap the bytes of the array elements
choose(choices[, out, mode])	Use an index array to construct a new array from a set of choices.
clip([min, max, out])	Return an array whose values are limited to [min, max].
compress(condition[, axis, out])	Return selected slices of this array along given axis.
conj()	Complex-conjugate all elements.
conjugate()	Return the complex conjugate, element-wise.
convert_units(new_unit)	Convert units of this array in-place.
copy([order])	Return a copy of the array.
cumprod([axis, dtype, out])	Return the cumulative product of the elements along the given axis.
cumsum([axis, dtype, out])	Return the cumulative sum of the elements along the given axis.
diagonal([offset, axis1, axis2])	Return specified diagonals.
dot(b[, out])	Dot product of two arrays.
dump(file)	Dump a pickle of the array to the specified file.
dumps()	Returns the pickle of the array as a string.
fill(value)	Fill the array with a scalar value.
flatten([order])	Return a copy of the array collapsed into one dimension.
getfield(dtype[, offset])	Returns a field of the given array as a certain type.
in_original_units()	Retun a copy of this array expressed in the units specified in the parameter file.
in_units(new_unit, **context_overrides)	Return a copy of this array expressed relative to an alternative unit.
item(*args)	Copy an element of an array to a standard Python scalar and return it.
itemset(*args)	Insert scalar into an array (scalar is cast to array’s dtype, if possible)
max([axis, out, keepdims, initial, where])	Return the maximum along a given axis.
mean([axis, dtype, out, keepdims])	Returns the average of the array elements along given axis.
min([axis, out, keepdims, initial, where])	Return the minimum along a given axis.
newbyteorder([new_order])	Return the array with the same data viewed with a different byte order.
nonzero()	Return the indices of the elements that are non-zero.
partition(kth[, axis, kind, order])	Rearranges the elements in the array in such a way that the value of the element in kth position is in the position it would be in a sorted array.
prod([axis, dtype, out, keepdims, initial, …])	Return the product of the array elements over the given axis
ptp([axis, out, keepdims])	Peak to peak (maximum - minimum) value along a given axis.
put(indices, values[, mode])	Set a.flat[n] = values[n] for all n in indices.
ravel([order])	Return a flattened array.
repeat(repeats[, axis])	Repeat elements of an array.
reshape(shape[, order])	Returns an array containing the same data with a new shape.
resize(new_shape[, refcheck])	Change shape and size of array in-place.
round([decimals, out])	Return a with each element rounded to the given number of decimals.
searchsorted(v[, side, sorter])	Find indices where elements of v should be inserted in a to maintain order.
set_default_units([quiet])	Set the units for this array by performing dimensional analysis on the default dimensions for the array.
set_units_like(new_unit)	Set the units for this array by performing dimensional analysis on the supplied unit and referring to the units of the original file
setfield(val, dtype[, offset])	Put a value into a specified place in a field defined by a data-type.
setflags([write, align, uic])	Set array flags WRITEABLE, ALIGNED, (WRITEBACKIFCOPY and UPDATEIFCOPY), respectively.
sort([axis, kind, order])	Sort an array in-place.
squeeze([axis])	Remove single-dimensional entries from the shape of a.
std([axis, dtype, out, ddof, keepdims])	Returns the standard deviation of the array elements along given axis.
sum([axis, dtype, out, keepdims, initial, where])	Return the sum of the array elements over the given axis.
swapaxes(axis1, axis2)	Return a view of the array with axis1 and axis2 interchanged.
take(indices[, axis, out, mode])	Return an array formed from the elements of a at the given indices.
tobytes([order])	Construct Python bytes containing the raw data bytes in the array.
tofile(fid[, sep, format])	Write array to a file as text or binary (default).
tolist()	Return the array as an a.ndim-levels deep nested list of Python scalars.
tostring([order])	Construct Python bytes containing the raw data bytes in the array.
trace([offset, axis1, axis2, dtype, out])	Return the sum along diagonals of the array.
transpose(*axes)	Returns a view of the array with axes transposed.
var([axis, dtype, out, ddof, keepdims])	Returns the variance of the array elements, along given axis.
view([dtype, type])	New view of array with the same data.
write(**kwargs)	Write this array to disk according to the standard method associated with its base file.

abs
conversion_context
mean_by_mass
ufunc_rule

property ancestor

Provides the basemost SimArray that an IndexedSimArray is based on.

cumsum(axis=None, dtype=None, out=None)

Return the cumulative sum of the elements along the given axis.

Refer to numpy.cumsum for full documentation.

See also

numpy.cumsum

equivalent function

prod(axis=None, dtype=None, out=None, keepdims=False, initial=1, where=True)

Return the product of the array elements over the given axis

Refer to numpy.prod for full documentation.

See also

numpy.prod

equivalent function

sum(axis=None, dtype=None, out=None, keepdims=False, initial=0, where=True)

Return the sum of the array elements over the given axis.

Refer to numpy.sum for full documentation.

See also

numpy.sum

equivalent function

mean(axis=None, dtype=None, out=None, keepdims=False)

Returns the average of the array elements along given axis.

Refer to numpy.mean for full documentation.

See also

numpy.mean

equivalent function

max(axis=None, out=None, keepdims=False, initial=<no value>, where=True)

Return the maximum along a given axis.

Refer to numpy.amax for full documentation.

See also

numpy.amax

equivalent function

min(axis=None, out=None, keepdims=False, initial=<no value>, where=True)

Return the minimum along a given axis.

Refer to numpy.amin for full documentation.

See also

numpy.amin

equivalent function

ptp(axis=None, out=None, keepdims=False)

Peak to peak (maximum - minimum) value along a given axis.

Refer to numpy.ptp for full documentation.

See also

numpy.ptp

equivalent function

std(axis=None, dtype=None, out=None, ddof=0, keepdims=False)

Returns the standard deviation of the array elements along given axis.

Refer to numpy.std for full documentation.

See also

numpy.std

equivalent function

var(axis=None, dtype=None, out=None, ddof=0, keepdims=False)

Returns the variance of the array elements, along given axis.

Refer to numpy.var for full documentation.

See also

numpy.var

equivalent function

set_units_like(new_unit)

Set the units for this array by performing dimensional analysis on the supplied unit and referring to the units of the original file

set_default_units(quiet=False)

Set the units for this array by performing dimensional analysis on the default dimensions for the array.

in_original_units()

Retun a copy of this array expressed in the units specified in the parameter file.

in_units(new_unit, **context_overrides)

Return a copy of this array expressed relative to an alternative unit.

convert_units(new_unit)

Convert units of this array in-place. Note that if this is a sub-view, the entire base array will be converted.

write(**kwargs)

Write this array to disk according to the standard method associated with its base file. This is equivalent to calling

>>> sim.gas.write_array('array')

in the case of writing out the array ‘array’ for the gas particle family. See the description of pynbody.snapshot.SimSnap.write_array() for options.

pynbody.array.wrapper_function(self, value, /)

Return self>=value.

exception pynbody.array.RemoteKeyboardInterrupt

pynbody.array.f(value, /)

Return self>=value.

pynbody.array.shared_array_remote(fn)

A decorator for functions returning a new function that is suitable for use remotely. Inputs to and outputs from the function can be transferred efficiently if they are backed onto shared memory. Ownership of any shared memory returned by the function is transferred.

pynbody.array.remote_map(pool, fn, *iterables)

A replacement for python’s in-built map function, sending out tasks to the pool and performing the magic required to transport shared memory arrays correctly. The function fn must be wrapped with the shared_array_remote decorator to interface correctly with this magic.

pynbody.array.exit_cleanup()

Clean up any shared memory that has not yet been freed. In theory this should not be required, but it is here as a safety net.

simdict

This submodule defines an augmented dictionary class (SimDict) for SimSnap properties where entries need to be managed e.g. for defining default entries, or for ensuring consistency between equivalent properties like redshift and scalefactor.

By default, a SimDict automatically converts between redshift (‘z’) and scalefactor (‘a’) and implements default entries for cosmological values listed in the [default-cosmology] section of the pynbody configuration files.

Adding further properties

To add further properties use the SimDict.getter and SimDict.setter decorators. For instance, to add a property ‘X_copy’ which just reflects the value of the property ‘X’, you would use the following code:

@SimDict.getter
def X_copy(d) :
    return d['X']

@SimDict.setter
def X_copy(d, value) :
    d['X'] = value

class pynbody.simdict.SimDict() -> new empty dictionary dict(mapping) -> new dictionary initialized from a mapping object's (key, value) pairs dict(iterable) -> new dictionary initialized as if via: d = {} for k, v in iterable: d[k] = v dict(**kwargs) -> new dictionary initialized with the name=value pairs in the keyword argument list. For example: dict(one=1, two=2)

Methods

clear()
copy()
fromkeys(iterable[, value])	Create a new dictionary with keys from iterable and values set to value.
get(key[, default])	Return the value for key if key is in the dictionary, else default.
getter(f)	Define a getter function for all SimDicts
items()
keys()
pop(k[,d])	If key is not found, d is returned if given, otherwise KeyError is raised
popitem(/)	Remove and return a (key, value) pair as a 2-tuple.
setdefault(key[, default])	Insert key with a value of default if key is not in the dictionary.
setter(f)	Define a setter function for all SimDicts
update([E, ]**F)	If E is present and has a .keys() method, then does: for k in E: D[k] = E[k] If E is present and lacks a .keys() method, then does: for k, v in E: D[k] = v In either case, this is followed by: for k in F: D[k] = F[k]
values()

static getter(f)

Define a getter function for all SimDicts

static setter(f)

Define a setter function for all SimDicts

family

This module defines the Family class which represents families of particles (e.g. dm, gas, star particles). New Family objects are automatically registered so that snapshots can use them in the normal syntax (snap.dm, snap.star, etc).

In practice the easiest way to make use of the flexibility this module provides is through adding more families of particles in your config.ini.

pynbody.family.family_names(with_aliases=False)

Returns a list of the names of all particle families. If with_aliases is True, include aliases in the list.

pynbody.family.get_family(name, create=False)

Returns a family corresponding to the specified string. If the family does not exist and create is False, raises ValueError. If the family does not exist and create is True, an appropriate object is instantiated, registered and returned.

snapshot

This module implements the SimSnap class which manages and stores snapshot data. It also implements the SubSnap class (and relatives) which represent different views of an existing SimSnap.

class pynbody.snapshot.SimSnap

The class for managing simulation snapshots.

For most purposes, SimSnaps should be initialized through load() or new().

For a basic tutorial explaining how to load a file as a SimSnap see tutorials/data_access.

Getting arrays or subsnaps

Once a SimSnap object f is instantiated, it can be used in various ways. The most common operation is to access something with the code f[x]. Depending on the type of x, various behaviours result:

• If x is a string, the array named by x is returned. If no such array exists, the framework attempts to load or derive an array of that name (in that order). If this is unsuccessful, a KeyError is raised.

• If x is a python slice (e.g. f[5:100:3]) or an array of integers (e.g. f[[1,5,100,200]]) a subsnap containing only the mentioned particles is returned.

  See tutorials/data_access for more information.

• If x is a numpy array of booleans, it is interpreted as a mask and a subsnap containing only those particles for which x[i] is True. This means that f[condition] is a shortcut for f[np.where(condition)].

• If x is a pynbody.filt.Filter object, a subsnap containing only the particles which pass the filter condition is returned.

  See tutorials/data_access for more information.

• If x is a pynbody.family.Family object, a subsnap containing only the particles in that family is returned. In practice for most code it is more convenient to write e.g. f.dm in place of the equivalent syntax f[pynbody.family.dm].

Getting metadata

The property filename gives the filename of a snapshot.

There is also a properties dictionary which contains further metadata about the snapshot. See Subsnaps.

Attributes

ancestor

The original SimSnap from which this view is derived (potentially self)

filename

Methods

all_keys()	Returns a list of all arrays that can be either lazy-evaluated or lazy loaded from an auxiliary file.
bridge(other)	Tries to construct a bridge function between this SimSnap and another one.
conversion_context()	Return a dictionary containing a (scalefactor) and h (Hubble constant in canonical units) for this snapshot, ready for passing into unit conversion functions.
derivable_keys()	Returns a list of arrays which can be lazy-evaluated.
derived_quantity(fn)
families()	Return the particle families which have representitives in this SimSnap.
family_keys([fam])	Return list of arrays which are not accessible from this view, but can be accessed from family-specific sub-views.
get(key[, alternative])	Standard python get method, returns self[key] if key in self else alternative
get_index_list(relative_to[, of_particles])	Get a list specifying the index of the particles in this view relative to the ancestor relative_to, such that relative_to[get_index_list(relative_to)]==self.
halos(*args, **kwargs)	Tries to instantiate a halo catalogue object for the given snapshot, using the first available method (as defined in the configuration files).
has_family_key(name)	Returns True if the array name is accessible (in memory) for at least one family
has_key(name)	Returns True if the array name is accessible (in memory)
infer_original_units(dimensions)	Given a unit (or string) dimensions, returns a unit with the same physical dimensions which is in the unit schema of the current file.
intersect(other[, op])	Returns the set intersection of this simulation view with another view of the same simulation
is_ancestor(other)	Returns true if other is a subview of self
is_derived_array(name[, fam])	Returns True if the array or family array of given name is auto-derived (and therefore read-only).
is_descendant(other)	Returns true if self is a subview of other
items()	Returns a list of tuples describing the array names and their contents in memory
keys()	Return the directly accessible array names (in memory)
load_copy()	Tries to load a copy of this snapshot, using partial loading to select only a subset of particles corresponding to a given SubSnap
loadable_keys([fam])	Returns a list of arrays which can be lazy-loaded from an auxiliary file.
mean_by_mass(name)	Calculate the mean by mass of the specified array.
original_units()	Converts all arrays’units to be consistent with the units of the original file.
physical_units([distance, velocity, mass, …])	Converts all array’s units to be consistent with the distance, velocity, mass basis units specified.
rotate_x(angle)	Rotates the snapshot about the current x-axis by ‘angle’ degrees.
rotate_y(angle)	Rotates the snapshot about the current y-axis by ‘angle’ degrees.
rotate_z(angle)	Rotates the snapshot about the current z-axis by ‘angle’ degrees.
set_units_system([velocity, distance, mass, …])	Set the unit system for the snapshot by specifying any or all of velocity, distance, mass and temperature units.
setdiff(other)	Returns the set difference of this simulation view with another view of the same simulation
transform(matrix)
union(other)	Returns the set union of this simulation view with another view of the same simulation
unlink_array(name)	If the named array is auto-derived, this destroys the link so that the array becomes editable but no longer auto-updates.
values()	Returns a list of the actual arrays in memory
wrap([boxsize, convention])	Wraps the positions of the particles in the box to lie between [-boxsize/2, boxsize/2].
write([fmt, filename])
write_array(array_name[, fam, overwrite])	Write out the array with the specified name.

decorator
iteritems
iterkeys
itervalues
stable_derived_quantity

families()

Return the particle families which have representitives in this SimSnap. The families are ordered by their appearance in the snapshot.

keys()

Return the directly accessible array names (in memory)

has_key(name)

Returns True if the array name is accessible (in memory)

values()

Returns a list of the actual arrays in memory

items()

Returns a list of tuples describing the array names and their contents in memory

get(key, alternative=None)

Standard python get method, returns self[key] if key in self else alternative

has_family_key(name)

Returns True if the array name is accessible (in memory) for at least one family

loadable_keys(fam=None)

Returns a list of arrays which can be lazy-loaded from an auxiliary file.

derivable_keys()

Returns a list of arrays which can be lazy-evaluated.

all_keys()

Returns a list of all arrays that can be either lazy-evaluated or lazy loaded from an auxiliary file.

family_keys(fam=None)

Return list of arrays which are not accessible from this view, but can be accessed from family-specific sub-views.

If fam is not None, only those keys applying to the specific family will be returned (equivalent to self.fam.keys()).

is_ancestor(other)

Returns true if other is a subview of self

is_descendant(other)

Returns true if self is a subview of other

property ancestor

The original SimSnap from which this view is derived (potentially self)

get_index_list(relative_to, of_particles=None)

Get a list specifying the index of the particles in this view relative to the ancestor relative_to, such that relative_to[get_index_list(relative_to)]==self.

intersect(other, op=<function intersect1d>)

Returns the set intersection of this simulation view with another view of the same simulation

union(other)

Returns the set union of this simulation view with another view of the same simulation

setdiff(other)

Returns the set difference of this simulation view with another view of the same simulation

conversion_context()

Return a dictionary containing a (scalefactor) and h (Hubble constant in canonical units) for this snapshot, ready for passing into unit conversion functions.

set_units_system(velocity=None, distance=None, mass=None, temperature=None)

Set the unit system for the snapshot by specifying any or all of velocity, distance, mass and temperature units. The units can be given as strings or as pynbody Unit objects.

If any of the units are not specified and a previous file_units_system does not exist, the defaults are used.

original_units()

Converts all arrays’units to be consistent with the units of the original file.

physical_units(distance='kpc', velocity='km s^-1', mass='Msol', persistent=True)

Converts all array’s units to be consistent with the distance, velocity, mass basis units specified.

Base units can be specified using keywords.

Optional Keywords:

distance: string (default = ‘kpc’)

velocity: string (default = ‘km s^-1’)

mass: string (default = ‘Msol’)

persistent: boolean (default = True); apply units change to future lazy-loaded arrays if True

infer_original_units(dimensions)

Given a unit (or string) dimensions, returns a unit with the same physical dimensions which is in the unit schema of the current file.

halos(*args, **kwargs)

Tries to instantiate a halo catalogue object for the given snapshot, using the first available method (as defined in the configuration files).

bridge(other)

Tries to construct a bridge function between this SimSnap and another one.

This function calls pynbody.bridge.bridge_factory(). For more information see bridge-tutorial, or the reference documentation for pynbody.bridge.

load_copy()

Tries to load a copy of this snapshot, using partial loading to select only a subset of particles corresponding to a given SubSnap

rotate_x(angle)

Rotates the snapshot about the current x-axis by ‘angle’ degrees.

rotate_y(angle)

Rotates the snapshot about the current y-axis by ‘angle’ degrees.

rotate_z(angle)

Rotates the snapshot about the current z-axis by ‘angle’ degrees.

wrap(boxsize=None, convention='center')

Wraps the positions of the particles in the box to lie between [-boxsize/2, boxsize/2].

If no boxsize is specified, self.properties[“boxsize”] is used.

write_array(array_name, fam=None, overwrite=False, **kwargs)

Write out the array with the specified name.

Some of the functionality is available via the pynbody.array.SimArray.write() method, which calls the present function with appropriate arguments.

Input

array_name - the name of the array to write

Optional Keywords

fam (None) - Write out only one family; or provide a list to

write out a set of families.

is_derived_array(name, fam=None)

Returns True if the array or family array of given name is auto-derived (and therefore read-only).

unlink_array(name)

If the named array is auto-derived, this destroys the link so that the array becomes editable but no longer auto-updates.

mean_by_mass(name)

Calculate the mean by mass of the specified array.

class pynbody.snapshot.SubSnap(base, _slice)

Represent a sub-view of a SimSnap, initialized by specifying a slice. Arrays accessed through __getitem__ are automatically sub-viewed using the given slice.

Attributes

ancestor

The original SimSnap from which this view is derived (potentially self)

filename

Methods

all_keys()	Returns a list of all arrays that can be either lazy-evaluated or lazy loaded from an auxiliary file.
bridge(other)	Tries to construct a bridge function between this SimSnap and another one.
conversion_context()	Return a dictionary containing a (scalefactor) and h (Hubble constant in canonical units) for this snapshot, ready for passing into unit conversion functions.
derivable_keys()	Returns a list of arrays which can be lazy-evaluated.
derived_quantity(fn)
families()	Return the particle families which have representitives in this SimSnap.
family_keys([fam])	Return list of arrays which are not accessible from this view, but can be accessed from family-specific sub-views.
get(key[, alternative])	Standard python get method, returns self[key] if key in self else alternative
get_index_list(relative_to[, of_particles])	Get a list specifying the index of the particles in this view relative to the ancestor relative_to, such that relative_to[get_index_list(relative_to)]==self.
halos(*args, **kwargs)	Tries to instantiate a halo catalogue object for the given snapshot, using the first available method (as defined in the configuration files).
has_family_key(name)	Returns True if the array name is accessible (in memory) for at least one family
has_key(name)	Returns True if the array name is accessible (in memory)
infer_original_units(*args)	Return the units on disk for a quantity with the specified dimensions
intersect(other[, op])	Returns the set intersection of this simulation view with another view of the same simulation
is_ancestor(other)	Returns true if other is a subview of self
is_derived_array(v[, fam])	Returns True if the array or family array of given name is auto-derived (and therefore read-only).
is_descendant(other)	Returns true if self is a subview of other
items()	Returns a list of tuples describing the array names and their contents in memory
keys()	Return the directly accessible array names (in memory)
load_copy()	Tries to load a copy of this snapshot, using partial loading to select only a subset of particles corresponding to a given SubSnap
loadable_keys([fam])	Returns a list of arrays which can be lazy-loaded from an auxiliary file.
mean_by_mass(name)	Calculate the mean by mass of the specified array.
original_units()	Converts all arrays’units to be consistent with the units of the original file.
physical_units(*args, **kwargs)	Converts all array’s units to be consistent with the distance, velocity, mass basis units specified.
rotate_x(angle)	Rotates the snapshot about the current x-axis by ‘angle’ degrees.
rotate_y(angle)	Rotates the snapshot about the current y-axis by ‘angle’ degrees.
rotate_z(angle)	Rotates the snapshot about the current z-axis by ‘angle’ degrees.
set_units_system([velocity, distance, mass, …])	Set the unit system for the snapshot by specifying any or all of velocity, distance, mass and temperature units.
setdiff(other)	Returns the set difference of this simulation view with another view of the same simulation
transform(matrix)
union(other)	Returns the set union of this simulation view with another view of the same simulation
unlink_array(name)	If the named array is auto-derived, this destroys the link so that the array becomes editable but no longer auto-updates.
values()	Returns a list of the actual arrays in memory
wrap([boxsize, convention])	Wraps the positions of the particles in the box to lie between [-boxsize/2, boxsize/2].
write([fmt, filename])
write_array(array_name[, fam])	Write out the array with the specified name.

decorator
iteritems
iterkeys
itervalues
stable_derived_quantity

keys()

Return the directly accessible array names (in memory)

loadable_keys(fam=None)

Returns a list of arrays which can be lazy-loaded from an auxiliary file.

derivable_keys()

Returns a list of arrays which can be lazy-evaluated.

infer_original_units(*args)

Return the units on disk for a quantity with the specified dimensions

write_array(array_name, fam=None, **kwargs)

Write out the array with the specified name.

Some of the functionality is available via the pynbody.array.SimArray.write() method, which calls the present function with appropriate arguments.

Input

array_name - the name of the array to write

Optional Keywords

fam (None) - Write out only one family; or provide a list to

write out a set of families.

family_keys(fam=None)

Return list of arrays which are not accessible from this view, but can be accessed from family-specific sub-views.

If fam is not None, only those keys applying to the specific family will be returned (equivalent to self.fam.keys()).

physical_units(*args, **kwargs)

Converts all array’s units to be consistent with the distance, velocity, mass basis units specified.

Base units can be specified using keywords.

Optional Keywords:

distance: string (default = ‘kpc’)

velocity: string (default = ‘km s^-1’)

mass: string (default = ‘Msol’)

persistent: boolean (default = True); apply units change to future lazy-loaded arrays if True

is_derived_array(v, fam=None)

Returns True if the array or family array of given name is auto-derived (and therefore read-only).

unlink_array(name)

If the named array is auto-derived, this destroys the link so that the array becomes editable but no longer auto-updates.

get_index_list(relative_to, of_particles=None)

Get a list specifying the index of the particles in this view relative to the ancestor relative_to, such that relative_to[get_index_list(relative_to)]==self.

class pynbody.snapshot.IndexedSubSnap(base, index_array=None, iord_array=None)

Represents a subset of the simulation particles according to an index array.

Parameters

baseSimSnap object

The base snapshot

index_arrayinteger array or None

The indices of the elements that define the sub snapshot. Set to None to use iord-based instead.

iord_arrayinteger array or None

The iord of the elements that define the sub snapshot. Set to None to use index-based instead. This may be computationally expensive. See note below.

Notes

index_array and iord_array arguments are mutually exclusive. In the case of iord_array, an sorting operation is required that may take a significant time and require O(N) memory.

Attributes

ancestor

The original SimSnap from which this view is derived (potentially self)

filename

Methods

all_keys()	Returns a list of all arrays that can be either lazy-evaluated or lazy loaded from an auxiliary file.
bridge(other)	Tries to construct a bridge function between this SimSnap and another one.
conversion_context()	Return a dictionary containing a (scalefactor) and h (Hubble constant in canonical units) for this snapshot, ready for passing into unit conversion functions.
derivable_keys()	Returns a list of arrays which can be lazy-evaluated.
derived_quantity(fn)
families()	Return the particle families which have representitives in this SimSnap.
family_keys([fam])	Return list of arrays which are not accessible from this view, but can be accessed from family-specific sub-views.
get(key[, alternative])	Standard python get method, returns self[key] if key in self else alternative
get_index_list(relative_to[, of_particles])	Get a list specifying the index of the particles in this view relative to the ancestor relative_to, such that relative_to[get_index_list(relative_to)]==self.
halos(*args, **kwargs)	Tries to instantiate a halo catalogue object for the given snapshot, using the first available method (as defined in the configuration files).
has_family_key(name)	Returns True if the array name is accessible (in memory) for at least one family
has_key(name)	Returns True if the array name is accessible (in memory)
infer_original_units(*args)	Return the units on disk for a quantity with the specified dimensions
intersect(other[, op])	Returns the set intersection of this simulation view with another view of the same simulation
is_ancestor(other)	Returns true if other is a subview of self
is_derived_array(v[, fam])	Returns True if the array or family array of given name is auto-derived (and therefore read-only).
is_descendant(other)	Returns true if self is a subview of other
items()	Returns a list of tuples describing the array names and their contents in memory
keys()	Return the directly accessible array names (in memory)
load_copy()	Tries to load a copy of this snapshot, using partial loading to select only a subset of particles corresponding to a given SubSnap
loadable_keys([fam])	Returns a list of arrays which can be lazy-loaded from an auxiliary file.
mean_by_mass(name)	Calculate the mean by mass of the specified array.
original_units()	Converts all arrays’units to be consistent with the units of the original file.
physical_units(*args, **kwargs)	Converts all array’s units to be consistent with the distance, velocity, mass basis units specified.
rotate_x(angle)	Rotates the snapshot about the current x-axis by ‘angle’ degrees.
rotate_y(angle)	Rotates the snapshot about the current y-axis by ‘angle’ degrees.
rotate_z(angle)	Rotates the snapshot about the current z-axis by ‘angle’ degrees.
set_units_system([velocity, distance, mass, …])	Set the unit system for the snapshot by specifying any or all of velocity, distance, mass and temperature units.
setdiff(other)	Returns the set difference of this simulation view with another view of the same simulation
transform(matrix)
union(other)	Returns the set union of this simulation view with another view of the same simulation
unlink_array(name)	If the named array is auto-derived, this destroys the link so that the array becomes editable but no longer auto-updates.
values()	Returns a list of the actual arrays in memory
wrap([boxsize, convention])	Wraps the positions of the particles in the box to lie between [-boxsize/2, boxsize/2].
write([fmt, filename])
write_array(array_name[, fam])	Write out the array with the specified name.

decorator
iteritems
iterkeys
itervalues
stable_derived_quantity

class pynbody.snapshot.FamilySubSnap(base, fam)

Represents a one-family portion of a parent snap object

Attributes

ancestor

The original SimSnap from which this view is derived (potentially self)

filename

Methods

all_keys()	Returns a list of all arrays that can be either lazy-evaluated or lazy loaded from an auxiliary file.
bridge(other)	Tries to construct a bridge function between this SimSnap and another one.
conversion_context()	Return a dictionary containing a (scalefactor) and h (Hubble constant in canonical units) for this snapshot, ready for passing into unit conversion functions.
derivable_keys()	Returns a list of arrays which can be lazy-evaluated.
derived_quantity(fn)
families()	Return the particle families which have representitives in this SimSnap.
family_keys([fam])	Return list of arrays which are not accessible from this view, but can be accessed from family-specific sub-views.
get(key[, alternative])	Standard python get method, returns self[key] if key in self else alternative
get_index_list(relative_to[, of_particles])	Get a list specifying the index of the particles in this view relative to the ancestor relative_to, such that relative_to[get_index_list(relative_to)]==self.
halos(*args, **kwargs)	Tries to instantiate a halo catalogue object for the given snapshot, using the first available method (as defined in the configuration files).
has_family_key(name)	Returns True if the array name is accessible (in memory) for at least one family
has_key(name)	Returns True if the array name is accessible (in memory)
infer_original_units(*args)	Return the units on disk for a quantity with the specified dimensions
intersect(other[, op])	Returns the set intersection of this simulation view with another view of the same simulation
is_ancestor(other)	Returns true if other is a subview of self
is_derived_array(v[, fam])	Returns True if the array or family array of given name is auto-derived (and therefore read-only).
is_descendant(other)	Returns true if self is a subview of other
items()	Returns a list of tuples describing the array names and their contents in memory
keys()	Return the directly accessible array names (in memory)
load_copy()	Tries to load a copy of this snapshot, using partial loading to select only a subset of particles corresponding to a given SubSnap
loadable_keys([fam])	Returns a list of arrays which can be lazy-loaded from an auxiliary file.
mean_by_mass(name)	Calculate the mean by mass of the specified array.
original_units()	Converts all arrays’units to be consistent with the units of the original file.
physical_units(*args, **kwargs)	Converts all array’s units to be consistent with the distance, velocity, mass basis units specified.
rotate_x(angle)	Rotates the snapshot about the current x-axis by ‘angle’ degrees.
rotate_y(angle)	Rotates the snapshot about the current y-axis by ‘angle’ degrees.
rotate_z(angle)	Rotates the snapshot about the current z-axis by ‘angle’ degrees.
set_units_system([velocity, distance, mass, …])	Set the unit system for the snapshot by specifying any or all of velocity, distance, mass and temperature units.
setdiff(other)	Returns the set difference of this simulation view with another view of the same simulation
transform(matrix)
union(other)	Returns the set union of this simulation view with another view of the same simulation
unlink_array(name)	If the named array is auto-derived, this destroys the link so that the array becomes editable but no longer auto-updates.
values()	Returns a list of the actual arrays in memory
wrap([boxsize, convention])	Wraps the positions of the particles in the box to lie between [-boxsize/2, boxsize/2].
write([fmt, filename])
write_array(array_name[, fam])	Write out the array with the specified name.

decorator
iteritems
iterkeys
itervalues
stable_derived_quantity

keys()

Return the directly accessible array names (in memory)

family_keys(fam=None)

Return list of arrays which are not accessible from this view, but can be accessed from family-specific sub-views.

If fam is not None, only those keys applying to the specific family will be returned (equivalent to self.fam.keys()).

pynbody.snapshot.load(filename, *args, **kwargs)

Loads a file using the appropriate class, returning a SimSnap instance.

pynbody.snapshot.new(n_particles=0, order=None, **families)

Create a blank SimSnap, with the specified number of particles.

Position, velocity and mass arrays are created and filled with zeros.

By default all particles are taken to be dark matter. To specify otherwise, pass in keyword arguments specifying the number of particles for each family, e.g.

f = new(dm=50, star=25, gas=25)

The order in which the different families appear in the snapshot is unspecified unless you add an ‘order’ argument:

f = new(dm=50, star=25, gas=25, order=’star,gas,dm’)

guarantees the stars, then gas, then dark matter particles appear in sequence.

units

The pynbody units module consists of a set of classes for tracking units.

It relates closely to the array module, which defines an extension to numpy arrays which carries unit information.

Making units

Units are generated and used at various points through the pynbody framework. Quite often the functions where users interact with units simply accept strings.

You can also make units yourself in two ways. Either you can create a string, and instantiate a Unit like this:

>>> units.Unit("Msol kpc**-3")
>>> units.Unit("2.1e12 m_p cm**-2/3")

Or you can do it within python, using the named Unit objects

>>> units.Msol * units.kpc**-3
>>> 2.1e12 * units.m_p * units.cm**(-2,3)

In the last example, either a tuple describing a fraction or a Fraction instance (from the standard python module fractions) is acceptable.

Getting conversion ratios

To convert one unit to another, use the ratio member function:

>>> units.Msol.ratio(units.kg)
1.99e30
>>> (units.Msol / units.kpc**3).ratio(units.m_p/units.cm**3)
4.04e-8

If the units cannot be converted, a UnitsException is raised:

>>> units.Msol.ratio(units.kpc)
UnitsException

Specifying numerical values

Sometimes it’s necessary to specify a numerical value in the course of a conversion. For instance, consider a comoving distance; this can be specified in pynbody units as follows:

>>> comoving_kpc = units.kpc * units.a

where units.a represents the scalefactor. We can attempt to convert this to a physical distance as follows

>>> comoving_kpc.ratio(units.kpc)

but this fails, throwing a UnitsException. On the other hand we can specify a value for the scalefactor when we request the conversion

>>> comoving_kpc.ratio(units.kpc, a=0.5)
0.5

and the conversion completes with the expected result. The units module also defines units.h for the dimensionless hubble constant, which can be used similarly. By default, all conversions happening within a specific simulation context should pass in values for a and h as a matter of routine.

Any IrreducibleUnit (see below) can have a value specified in this way, but a and h are envisaged to be the most useful applications.

Defining new base units

The units module is fully extensible: you can define and name your own units which then integrate with all the standard functions.

litre = units.NamedUnit("litre",0.001*units.m**3)
gallon = units.NamedUnit("gallon",0.004546*units.m**3)
gallon.ratio(litre) # 4.546
(units.pc**3).ratio(litre) # 2.94e52

You can even define completely new dimensions.

V = units.IrreducibleUnit("V") # define a volt
C = units.NamedUnit("C", units.J/V) # define a coulomb
q = units.NamedUnit("q", 1.60217646e-19*C) # elementary charge
F = units.NamedUnit("F", C/V) # Farad
epsilon0 = 8.85418e-12 *F/units.m

>>> (q*V).ratio("eV")
1.000
>>> ((q**2)/(4*math.pi*epsilon0*units.m**2)).ratio("N")
2.31e-28

exception pynbody.units.UnitsException

class pynbody.units.UnitBase

Base class for units. To instantiate a unit, call the pynbody.units.Unit() factory function.

Methods

in_units(*a, **kw)	Alias for ratio
irrep()	Return a unit equivalent to this one (may be identical) but expressed in terms of the currently defined IrreducibleUnit instances.
ratio(other, **substitutions)	Get the conversion ratio between this Unit and another specified unit.

is_dimensionless
simplify

ratio(other, **substitutions)

Get the conversion ratio between this Unit and another specified unit.

Keyword arguments, if specified, give numerical substitutions for the named unit. This is most useful for specifying values for cosmological quantities like ‘a’ and ‘h’, but can also be used for any IrreducibleUnit.

>>> Unit("1 Mpc a").ratio("kpc", a=0.25)
250.0
>>> Unit("1 Mpc").ratio("Msol")
UnitsException: not convertible
>>> Unit("1 Mpc").ratio("Msol", kg=25.0, m=50.0)
3.1028701506345152e-08

in_units(*a, **kw)

Alias for ratio

irrep()

Return a unit equivalent to this one (may be identical) but expressed in terms of the currently defined IrreducibleUnit instances.

class pynbody.units.NoUnit

Methods

in_units(*a, **kw)	Alias for ratio
irrep()	Return a unit equivalent to this one (may be identical) but expressed in terms of the currently defined IrreducibleUnit instances.
ratio(other, **substitutions)	Get the conversion ratio between this Unit and another specified unit.

dimensional_project
is_dimensionless
latex
simplify

ratio(other, **substitutions)

Get the conversion ratio between this Unit and another specified unit.

Keyword arguments, if specified, give numerical substitutions for the named unit. This is most useful for specifying values for cosmological quantities like ‘a’ and ‘h’, but can also be used for any IrreducibleUnit.

>>> Unit("1 Mpc a").ratio("kpc", a=0.25)
250.0
>>> Unit("1 Mpc").ratio("Msol")
UnitsException: not convertible
>>> Unit("1 Mpc").ratio("Msol", kg=25.0, m=50.0)
3.1028701506345152e-08

irrep()

Return a unit equivalent to this one (may be identical) but expressed in terms of the currently defined IrreducibleUnit instances.

class pynbody.units.IrreducibleUnit(st)

Methods

in_units(*a, **kw)	Alias for ratio
irrep()	Return a unit equivalent to this one (may be identical) but expressed in terms of the currently defined IrreducibleUnit instances.
ratio(other, **substitutions)	Get the conversion ratio between this Unit and another specified unit.

is_dimensionless
latex
simplify

irrep()

Return a unit equivalent to this one (may be identical) but expressed in terms of the currently defined IrreducibleUnit instances.

class pynbody.units.NamedUnit(st, represents)

Methods

in_units(*a, **kw)	Alias for ratio
irrep()	Return a unit equivalent to this one (may be identical) but expressed in terms of the currently defined IrreducibleUnit instances.
ratio(other, **substitutions)	Get the conversion ratio between this Unit and another specified unit.

is_dimensionless
latex
simplify

irrep()

Return a unit equivalent to this one (may be identical) but expressed in terms of the currently defined IrreducibleUnit instances.

class pynbody.units.CompositeUnit(scale, bases, powers)

Methods

copy()	Create a copy which is ‘shallow’ in the sense that it references exactly the same underlying base units, but where the list of those units can be manipulated separately.
dimensional_project(basis_units)	Work out how to express the dimensions of this unit relative to the specified list of basis units.
dimensionless_constant(**substitutions)	If this unit is dimensionless, return its scalar quantity.
in_units(*a, **kw)	Alias for ratio
irrep()	Return a new unit which represents this unit expressed solely in terms of IrreducibleUnit bases.
is_dimensionless()	Returns true if this unit actually translates into a scalar quantity.
latex()	Returns a LaTeX representation of this unit.
ratio(other, **substitutions)	Get the conversion ratio between this Unit and another specified unit.

simplify

latex()

Returns a LaTeX representation of this unit.

Prefactors are converted into exponent notation. Named units by default are represented by the string ‘mathrm{unit_name}’, although this can be overriden in the pynbody configuration files or by setting unit_name._latex.

copy()

Create a copy which is ‘shallow’ in the sense that it references exactly the same underlying base units, but where the list of those units can be manipulated separately.

irrep()

Return a new unit which represents this unit expressed solely in terms of IrreducibleUnit bases.

is_dimensionless()

Returns true if this unit actually translates into a scalar quantity.

dimensionless_constant(**substitutions)

If this unit is dimensionless, return its scalar quantity.

Direct use of this function is not recommended. It is generally better to use the ratio function instead.

Provide keyword arguments to set values for named IrreducibleUnits – see the ratio function for more information.

dimensional_project(basis_units)

Work out how to express the dimensions of this unit relative to the specified list of basis units.

This is used by the framework when making inferences about sensible units to use in various situations.

For example, you can represent a length as an energy divided by a force:

>>> Unit("23 kpc").dimensional_project(["J", "N"])
array([1, -1], dtype=object)

However it’s not possible to represent a length by energy alone:

>>> Unit("23 kpc").dimensional_project(["J"])
UnitsException: Basis units do not span dimensions of specified unit

This function also doesn’t know what to do if the result is ambiguous:

>>> Unit("23 kpc").dimensional_project(["J", "N", "kpc"])
UnitsException: Basis units are not linearly independent

pynbody.units.Unit(s)

Class factory for units. Given a string s, creates a Unit object.

The string format is:

[<scale>] [<unit_name>][**<rational_power>] [[<unit_name>] … ]

for example:

“1.e30 kg”

“kpc**2”

“26.2 m s**-1”

pynbody.units.takes_arg_in_units(*args, **orig_kwargs)

Returns a decorator to create a function which auto-converts input to given units.

Usage:

@takes_arg_in_units((2, "Msol"), (1, "kpc"), ("blob", "erg"))
def my_function(arg0, arg1, arg2, blob=22) :
   print "Arg 2 is",arg2,"Msol"
   print "Arg 1 is",arg1,"kpc"
   print "blob is",blob,"ergs"

>>> My_function(22, "1.e30 kg", 23, blob="12 J")
Input 3 is 0.5 Msol
Input 2 is 23 kpc

pynbody.units.has_unit(obj)

Returns True if the specified object has a meaningful units attribute

pynbody.units.has_units(obj)

Returns True if the specified object has a meaningful units attribute

pynbody.units.is_unit(obj)

Returns True if the specified object represents a unit

pynbody.units.is_unit_like(obj)

Returns True if the specified object is itself a unit or otherwise exposes unit information

util

Various utility routines used internally by pynbody.

pynbody.util.open_(filename, *args)

Open a file, determining from the filename whether to use gzip decompression

pynbody.util.open_with_size(filename, *args)

Open a file for reading, returning also the (decompressed) file size

pynbody.util.eps_as_simarray(f, eps)

Convert th given eps to a SimArray with units of f[‘pos’] and dtype of f[‘mass’]

pynbody.util.get_eps(f)

The gravitational softening length is determined from (in order of preference): 1. the array f[‘eps’] 2. f.properties[‘eps’] (scalar or unit)

Return a SimArray with correct units and dtype (same dtype as ‘mass’ array)

pynbody.util.intersect_slices(s1, s2, array_length=None)

Given two python slices s1 and s2, return a new slice which will extract the data of an array d which is in both d[s1] and d[s2].

Note that it may not be possible to do this without information on the length of the array referred to, hence all slices with end-relative indexes are first converted into begin-relative indexes. This means that the slice returned may be specific to the length specified.

pynbody.util.relative_slice(s_relative_to, s)

Given a slice s, return a slice s_prime with the property that array[s_relative_to][s_prime] == array[s]. Clearly this will not be possible for arbitrarily chosen s_relative_to and s, but it should be possible for s=intersect_slices(s_relative_to, s_any) which is the use case envisioned here (and used by SubSim). This code currently does not work with end-relative (i.e. negative) start or stop positions.

pynbody.util.chained_slice(s1, s2)

Return a slice s3 with the property that ar[s1][s2] == ar[s3]

pynbody.util.index_before_slice(s, index)

Return an index array new_index with the property that, for a slice s (start, stop and step all positive), ar[s][index] == ar[new_index].

pynbody.util.concatenate_indexing(i1, i2)

Given either a numpy array or slice for both i1 and i2, return either a numpy array or slice i3 with the property that

ar[i3] == ar[i1][i2].

As a convenience, if i2 is None, i1 is returned

pynbody.util.indexing_length(sl_or_ar)

Given either an array or slice, return len(ar[sl_or_ar]) for any array ar which is large enough that the slice does not overrun it.

pynbody.util.arrays_are_same(a1, a2)

Returns True if a1 and a2 are numpy views pointing to the exact same underlying data; False otherwise.

pynbody.util.set_array_if_not_same(a_store, a_in, index=None)

This routine checks whether a_store and a_in ultimately point to the same buffer; if not, the contents of a_in are copied into a_store.

pynbody.util.index_of_first(array, find)

Returns the index to the first element in array which satisfies array[index]>=find. The array must be sorted in ascending order.

pynbody.util.equipartition(ar, nbins, vmin=None, vmax=None)

Given an array ar, return nbins+1 monotonically increasing bin edges such that the number of items in each bin is approximately equal.

pynbody.util.bisect(left, right, f, epsilon=None, eta=0, verbose=False, niter_max=200)

Finds the value x such that f(x)=0 for a monotonically increasing function f, using a binary search.

The search stops when either the bounding domain is smaller than epsilon (by default 10^-7 times the original region) OR a value f(x) is found such that |f(x)|<eta (by default eta=0, so this criterion is never satisfied).

pynbody.util.gauss_jordan(out)

A simple Gauss-Jordan matrix inverter. This is provided so that matrices of fractions can be inverted (numpy linalg converts everything to floats first.)

Don’t use on large matrices – it’s slow!

Based on public domain code by Jarno Elonen.

pynbody.util.rational_matrix_inv(matrix)

A simple replacement for numpy linalg matrix inverse which handles fractions exactly. Not suitable for large matrices!

pynbody.util.random_rotation_matrix()

Return a random rotation matrix (Haar measure for 3x3 case), using fast algorithm from Graphics Gems III

(http://tog.acm.org/resources/GraphicsGems/gemsiii/rand_rotation.c)

pynbody.util.cutgz(x)

Strip the .gz ending off a string

pynbody.util.gamma_inc(a, z, eps=3e-07)

Incomplete gamma function accepting complex z, based on algorithm given in numerical recipes (3rd ed)

pynbody.util.threadsafe_inline(*args, **kwargs)

When scipy.weave.inline is called, it may trigger a compile. We only want one compilation to be going on at once, otherwise nasty race conditions arise. This function wraps scipy.weave.inline to be thread-safe.