API Documentation

pytreegrav.frontend.Accel(pos, m, softening=None, G=1.0, theta=0.7, tree=None, return_tree=False, parallel=False, method='adaptive', quadrupole=False)

Gravitational acceleration calculation

Returns the gravitational acceleration for a set of particles with positions x and masses m, at the positions of those particles, using either brute force or tree-based methods depending on the number of particles.

Parameters
  • pos (array_like) – shape (N,3) array of particle positions

  • m (array_like) – shape (N,) array of particle masses

  • G (float, optional) – gravitational constant (default 1.0)

  • softening (None or array_like, optional) – shape (N,) array containing kernel support radii for gravitational softening - these give the radius of compact support of the M4 cubic spline mass distribution - set to 0 by default

  • theta (float, optional) – cell opening angle used to control force accuracy; smaller is slower (runtime ~ theta^-3) but more accurate. (default 0.7, gives ~1% accuracy)

  • parallel (bool, optional) – If True, will parallelize the force summation over all available cores. (default False)

  • tree (Octree, optional) – optional pre-generated Octree: this can contain any set of particles, not necessarily the target particles at pos (default None)

  • return_tree (bool, optional) – return the tree used for future use (default False)

  • method (str, optional) – Which summation method to use: ‘adaptive’, ‘tree’, or ‘bruteforce’ (default adaptive tries to pick the faster choice)

  • quadrupole (bool, optional) – Whether to use quadrupole moments in tree summation (default False)

Returns

g – shape (N,3) array of acceleration vectors at the particle positions

Return type

array_like

pytreegrav.frontend.AccelTarget(pos_target, pos_source, m_source, softening_target=None, softening_source=None, G=1.0, theta=0.7, tree=None, return_tree=False, parallel=False, method='adaptive', quadrupole=False)

Gravitational acceleration calculation for general N+M body case

Returns the gravitational acceleration for a set of M particles with positions x_source and masses m_source, at the positions of a set of N particles that need not be the same.

Parameters
  • pos_target (array_like) – shape (N,3) array of target particle positions where you want to know the acceleration

  • pos_source (array_like) – shape (M,3) array of source particle positions (positions of particles sourcing the gravitational field)

  • m_source (array_like) – shape (M,) array of source particle masses

  • softening_target (array_like or None, optional) – shape (N,) array of target particle softening radii - these give the radius of compact support of the M4 cubic spline mass distribution

  • softening_source (array_like or None, optional) – shape (M,) array of source particle radii - these give the radius of compact support of the M4 cubic spline mass distribution

  • G (float, optional) – gravitational constant (default 1.0)

  • theta (float, optional) – cell opening angle used to control force accuracy; smaller is slower (runtime ~ theta^-3) but more accurate. (default 0.7, gives ~1% accuracy)

  • parallel (bool, optional) – If True, will parallelize the force summation over all available cores. (default False)

  • tree (Octree, optional) – optional pre-generated Octree: this can contain any set of particles, not necessarily the target particles at pos (default None)

  • return_tree (bool, optional) – return the tree used for future use (default False)

  • method (str, optional) – Which summation method to use: ‘adaptive’, ‘tree’, or ‘bruteforce’ (default adaptive tries to pick the faster choice)

  • quadrupole (bool, optional) – Whether to use quadrupole moments in tree summation (default False)

Returns

phi – shape (N,3) array of accelerations at the target positions

Return type

array_like

pytreegrav.frontend.ColumnDensity(pos, m, radii, rays=None, tree=None, return_tree=False, parallel=False)

Ray-traced column density calculation.

Returns the column density from the position of each particle integrated to infinity, assuming the particles are represented by uniform spheres. Note that optical depth can be obtained by supplying “σ = opacity * mass” in place of mass, useful in situations where opacity is highly variable.

Parameters
  • pos (array_like) – shape (N,3) array of particle positions

  • m (array_like) – shape (N,) array of particle masses

  • radii (array_like) – shape (N,) array containing particle radii of the uniform spheres that we use to model the particles’ mass distribution

  • rays (optional) – Which ray directions to raytrace the columns. DEFAULT: The simple 6-ray grid. OPTION 2: Give a (N_rays,3) array of vectors specifying the directions, which will be automatically normalized. OPTION 3: Give an integer number N to generate a raygrid of N random directions.

  • parallel (bool, optional) – If True, will parallelize the column density over all available cores. (default False)

  • tree (Octree, optional) – optional pre-generated Octree: this can contain any set of particles, not necessarily the target particles at pos (default None)

  • return_tree (bool, optional) – return the tree used for future use (default False)

Returns

columns – shape (N,N_rays) float array of column densities from particle centers integrated along the rays

Return type

array_like

pytreegrav.frontend.ConstructTree(pos, m, softening=None, quadrupole=False, vel=None)

Builds a tree containing particle data, for subsequent potential/field evaluation

Parameters
  • pos (array_like) – shape (N,3) array of particle positions

  • m (array_like) – shape (N,) array of particle masses

  • softening (array_like or None, optional) – shape (N,) array of particle softening lengths - these give the radius of compact support of the M4 cubic spline mass distribution of each particle

  • quadrupole (bool, optional) – Whether to store quadrupole moments (default False)

  • vel (bool, optional) – Whether to store node velocities in the tree (default False)

Returns

tree – Octree instance built from particle data

Return type

octree

pytreegrav.frontend.DensityCorrFunc(pos, m, rbins=None, max_bin_size_ratio=100, theta=1.0, tree=None, return_tree=False, parallel=False, boxsize=0, weighted_binning=False)

Computes the average amount of mass in radial bin [r,r+dr] around a point, provided a set of radial bins.

Parameters
  • pos (array_like) – shape (N,3) array of particle positions

  • m (array_like) – shape (N,) array of particle masses

  • rbins (array_like or None, optional) – 1D array of radial bin edges - if None will use heuristics to determine sensible bins. Otherwise MUST BE LOGARITHMICALLY SPACED (default None)

  • max_bin_size_ratio (float, optional) – controls the accuracy of the binning - tree nodes are subdivided until their side length is at most this factor * the radial bin width (default 100)

  • theta (float, optional) – cell opening angle used to control accuracy; smaller is slower (runtime ~ theta^-3) but more accurate. (default 1.0)

  • parallel (bool, optional) – If True, will parallelize the correlation function computation over all available cores. (default False)

  • tree (Octree, optional) – optional pre-generated Octree: this can contain any set of particles, not necessarily the target particles at pos (default None)

  • return_tree (bool, optional) – if True will return the generated or used tree for future use (default False)

  • boxsize (float, optional) – finite periodic box size, if periodic boundary conditions are to be used (default 0)

  • weighted_binning (bool, optional) – (experimental) if True will distribute mass among radial bings with a weighted kernel (default False)

Returns

  • rbins (array_like) – array containing radial bin edges

  • mbins (array_like) – array containing mean mass in radial bins, averaged over all points

pytreegrav.frontend.Potential(pos, m, softening=None, G=1.0, theta=0.7, tree=None, return_tree=False, parallel=False, method='adaptive', quadrupole=False)

Gravitational potential calculation

Returns the gravitational potential for a set of particles with positions x and masses m, at the positions of those particles, using either brute force or tree-based methods depending on the number of particles.

Parameters
  • pos (array_like) – shape (N,3) array of particle positions

  • m (array_like) – shape (N,) array of particle masses

  • G (float, optional) – gravitational constant (default 1.0)

  • softening (None or array_like, optional) – shape (N,) array containing kernel support radii for gravitational softening - - these give the radius of compact support of the M4 cubic spline mass distribution - set to 0 by default

  • theta (float, optional) – cell opening angle used to control force accuracy; smaller is slower (runtime ~ theta^-3) but more accurate. (default 0.7, gives ~1% accuracy)

  • parallel (bool, optional) – If True, will parallelize the force summation over all available cores. (default False)

  • tree (Octree, optional) – optional pre-generated Octree: this can contain any set of particles, not necessarily the target particles at pos (default None)

  • return_tree (bool, optional) – return the tree used for future use (default False)

  • method (str, optional) – Which summation method to use: ‘adaptive’, ‘tree’, or ‘bruteforce’ (default adaptive tries to pick the faster choice)

  • quadrupole (bool, optional) – Whether to use quadrupole moments in tree summation (default False)

Returns

phi – shape (N,) array of potentials at the particle positions

Return type

array_like

pytreegrav.frontend.PotentialTarget(pos_target, pos_source, m_source, softening_target=None, softening_source=None, G=1.0, theta=0.7, tree=None, return_tree=False, parallel=False, method='adaptive', quadrupole=False)

Gravitational potential calculation for general N+M body case

Returns the gravitational potential for a set of M particles with positions x_source and masses m_source, at the positions of a set of N particles that need not be the same.

Parameters
  • pos_target (array_like) – shape (N,3) array of target particle positions where you want to know the potential

  • pos_source (array_like) – shape (M,3) array of source particle positions (positions of particles sourcing the gravitational field)

  • m_source (array_like) – shape (M,) array of source particle masses

  • softening_target (array_like or None, optional) – shape (N,) array of target particle softening radii - these give the radius of compact support of the M4 cubic spline mass distribution

  • softening_source (array_like or None, optional) – shape (M,) array of source particle radii - these give the radius of compact support of the M4 cubic spline mass distribution

  • G (float, optional) – gravitational constant (default 1.0)

  • theta (float, optional) – cell opening angle used to control force accuracy; smaller is slower (runtime ~ theta^-3) but more accurate. (default 0.7, gives ~1% accuracy)

  • parallel (bool, optional) – If True, will parallelize the force summation over all available cores. (default False)

  • tree (Octree, optional) – optional pre-generated Octree: this can contain any set of particles, not necessarily the target particles at pos (default None)

  • return_tree (bool, optional) – return the tree used for future use (default False)

  • method (str, optional) – Which summation method to use: ‘adaptive’, ‘tree’, or ‘bruteforce’ (default adaptive tries to pick the faster choice)

  • quadrupole (bool, optional) – Whether to use quadrupole moments in tree summation (default False)

Returns

phi – shape (N,) array of potentials at the target positions

Return type

array_like

pytreegrav.frontend.VelocityCorrFunc(pos, m, v, rbins=None, max_bin_size_ratio=100, theta=1.0, tree=None, return_tree=False, parallel=False, boxsize=0, weighted_binning=False)

Computes the weighted average product v(x).v(x+r), for a vector field v, in radial bins

Parameters
  • pos (array_like) – shape (N,3) array of particle positions

  • m (array_like) – shape (N,) array of particle masses

  • v (array_like) – shape (N,3) of vector quantity (e.g. velocity, magnetic field, etc)

  • rbins (array_like or None, optional) – 1D array of radial bin edges - if None will use heuristics to determine sensible bins. Otherwise MUST BE LOGARITHMICALLY SPACED (default None)

  • max_bin_size_ratio (float, optional) – controls the accuracy of the binning - tree nodes are subdivided until their side length is at most this factor * the radial bin width (default 100)

  • theta (float, optional) – cell opening angle used to control accuracy; smaller is slower (runtime ~ theta^-3) but more accurate. (default 1.0)

  • parallel (bool, optional) – If True, will parallelize the correlation function computation over all available cores. (default False)

  • tree (Octree, optional) – optional pre-generated Octree: this can contain any set of particles, not necessarily the target particles at pos (default None)

  • return_tree (bool, optional) – if True will return the generated or used tree for future use (default False)

  • boxsize (float, optional) – finite periodic box size, if periodic boundary conditions are to be used (default 0)

  • weighted_binning (bool, optional) – (experimental) if True will distribute mass among radial bings with a weighted kernel (default False)

Returns

  • rbins (array_like) – array containing radial bin edges

  • corr (array_like) – array containing correlation function values in radial bins

pytreegrav.frontend.VelocityStructFunc(pos, m, v, rbins=None, max_bin_size_ratio=100, theta=1.0, tree=None, return_tree=False, parallel=False, boxsize=0, weighted_binning=False)

Computes the structure function for a vector field: the average value of |v(x)-v(x+r)|^2, in radial bins for r

Parameters
  • pos (array_like) – shape (N,3) array of particle positions

  • m (array_like) – shape (N,) array of particle masses

  • v (array_like) – shape (N,3) of vector quantity (e.g. velocity, magnetic field, etc)

  • rbins (array_like or None, optional) – 1D array of radial bin edges - if None will use heuristics to determine sensible bins. Otherwise MUST BE LOGARITHMICALLY SPACED (default None)

  • max_bin_size_ratio (float, optional) – controls the accuracy of the binning - tree nodes are subdivided until their side length is at most this factor * the radial bin width (default 100)

  • theta (float, optional) – cell opening angle used to control accuracy; smaller is slower (runtime ~ theta^-3) but more accurate. (default 1.0)

  • parallel (bool, optional) – If True, will parallelize the correlation function computation over all available cores. (default False)

  • tree (Octree, optional) – optional pre-generated Octree: this can contain any set of particles, not necessarily the target particles at pos (default None)

  • return_tree (bool, optional) – if True will return the generated or used tree for future use (default False)

  • boxsize (float, optional) – finite periodic box size, if periodic boundary conditions are to be used (default 0)

  • weighted_binning (bool, optional) – (experimental) if True will distribute mass among radial bings with a weighted kernel (default False)

Returns

  • rbins (array_like) – array containing radial bin edges

  • Sv (array_like) – array containing structure function values