pytreegrav package
Submodules
pytreegrav.bruteforce module
-
pytreegrav.bruteforce.
AccelTarget_bruteforce
(x_target, softening_target, x_source, m_source, softening_source, G=1.0) Returns the gravitational acceleration at a set of target positions, due to a set of source particles.
Arguments: x_target – shape (N,3) array of positions where the field is to be evaluated softening_target – shape (N,) array of minimum softening lengths to be used x_source – shape (M,3) array of positions of gravitating particles m_source – shape (M,) array of particle masses softening_source – shape (M,) array of softening lengths
Optional arguments: G – gravitational constant (default 1.0)
Returns: shape (N,3) array of gravitational accelerations
-
pytreegrav.bruteforce.
AccelTarget_bruteforce_parallel
(x_target, softening_target, x_source, m_source, softening_source, G=1.0) Returns the gravitational acceleration at a set of target positions, due to a set of source particles.
Arguments: x_target – shape (N,3) array of positions where the field is to be evaluated softening_target – shape (N,) array of minimum softening lengths to be used x_source – shape (M,3) array of positions of gravitating particles m_source – shape (M,) array of particle masses softening_source – shape (M,) array of softening lengths
Optional arguments: G – gravitational constant (default 1.0)
Returns: shape (N,3) array of gravitational accelerations
-
pytreegrav.bruteforce.
Accel_bruteforce
(x, m, softening, G=1.0) Returns the exact mutually-interacting gravitational accelerations of a set of particles.
Arguments: x – shape (N,3) array of positions where the potential is to be evaluated m – shape (N,) array of particle masses softening – shape (N,) array of softening lengths
Optional arguments: G – gravitational constant (default 1.0)
Returns: shape (N,3) array of gravitational accelerations
-
pytreegrav.bruteforce.
Accel_bruteforce_parallel
(x, m, softening, G=1.0) Returns the exact mutually-interacting gravitational accelerations of a set of particles.
Arguments: x – shape (N,3) array of positions where the potential is to be evaluated m – shape (N,) array of particle masses softening – shape (N,) array of softening lengths
Optional arguments: G – gravitational constant (default 1.0)
Returns: shape (N,3) array of gravitational accelerations
-
pytreegrav.bruteforce.
PotentialTarget_bruteforce
(x_target, softening_target, x_source, m_source, softening_source, G=1.0) Returns the exact gravitational potential due to a set of particles, at a set of positions that need not be the same as the particle positions.
Arguments: x_target – shape (N,3) array of positions where the potential is to be evaluated softening_target – shape (N,) array of minimum softening lengths to be used x_source – shape (M,3) array of positions of gravitating particles m_source – shape (M,) array of particle masses softening_source – shape (M,) array of softening lengths
Optional arguments: G – gravitational constant (default 0.7)
Returns: shape (N,) array of potential values
-
pytreegrav.bruteforce.
PotentialTarget_bruteforce_parallel
(x_target, softening_target, x_source, m_source, softening_source, G=1.0) Returns the exact gravitational potential due to a set of particles, at a set of positions that need not be the same as the particle positions.
Arguments: x_target – shape (N,3) array of positions where the potential is to be evaluated softening_target – shape (N,) array of minimum softening lengths to be used x_source – shape (M,3) array of positions of gravitating particles m_source – shape (M,) array of particle masses softening_source – shape (M,) array of softening lengths
Optional arguments: G – gravitational constant (default 0.7)
Returns: shape (N,) array of potential values
-
pytreegrav.bruteforce.
Potential_bruteforce
(x, m, softening, G=1.0) Returns the exact mutually-interacting gravitational potential for a set of particles with positions x and masses m, evaluated by brute force.
Arguments: x – shape (N,3) array of particle positions m – shape (N,) array of particle masses softening – shape (N,) array containing kernel support radii for gravitational softening
Optional arguments: G – gravitational constant (default 1.0)
Returns: shape (N,) array containing potential values
-
pytreegrav.bruteforce.
Potential_bruteforce_parallel
(x, m, softening, G=1.0) Returns the exact mutually-interacting gravitational potential for a set of particles with positions x and masses m, evaluated by brute force.
Arguments: x – shape (N,3) array of particle positions m – shape (N,) array of particle masses softening – shape (N,) array containing kernel support radii for gravitational softening
Optional arguments: G – gravitational constant (default 1.0)
Returns: shape (N,) array containing potential values
pytreegrav.dynamic_tree module
-
pytreegrav.dynamic_tree.
ComputeMomentsDynamic
(tree, no, children)
-
pytreegrav.dynamic_tree.
SetupTreewalk
(tree, no, children)
pytreegrav.frontend module
-
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
-
pytreegrav.frontend.
valueTestMethod
(method)
pytreegrav.kernel module
-
pytreegrav.kernel.
ForceKernel
(r, h) Returns the quantity equivalent to (fraction of mass enclosed)/ r^3 for a cubic-spline mass distribution of compact support radius h. Used to calculate the softened gravitational force.
Arguments: r - radius h - softening
-
pytreegrav.kernel.
PotentialKernel
(r, h) Returns the equivalent of -1/r for a cubic-spline mass distribution of compact support radius h. Used to calculate the softened gravitational potential.
Arguments: r - radius h - softening
pytreegrav.octree module
-
pytreegrav.octree.
ComputeMoments
(tree, no, children)
-
pytreegrav.octree.
SetupTreewalk
(tree, no, children)
pytreegrav.treewalk module
-
pytreegrav.treewalk.
AccelTarget_tree
(pos_target, softening_target, tree, G=1.0, theta=0.7, quadrupole=False) Returns the gravitational acceleration at the specified points, given a tree containing the mass distribution Arguments: pos_target – shape (N,3) array of positions at which to evaluate the field softening_target – shape (N,) array of minimum softening lengths to be used in all accel computations tree – Octree instance containing the positions, masses, and softenings of the source particles Optional arguments: G – gravitational constant (default 1.0) theta – accuracy parameter, smaller is more accurate, larger is faster (default 0.7) Returns: shape (N,3) array of acceleration values at each point in pos_target
-
pytreegrav.treewalk.
AccelTarget_tree_parallel
(pos_target, softening_target, tree, G=1.0, theta=0.7, quadrupole=False) Returns the gravitational acceleration at the specified points, given a tree containing the mass distribution Arguments: pos_target – shape (N,3) array of positions at which to evaluate the field softening_target – shape (N,) array of minimum softening lengths to be used in all accel computations tree – Octree instance containing the positions, masses, and softenings of the source particles Optional arguments: G – gravitational constant (default 1.0) theta – accuracy parameter, smaller is more accurate, larger is faster (default 0.7) Returns: shape (N,3) array of acceleration values at each point in pos_target
-
pytreegrav.treewalk.
AccelWalk
(pos, tree, softening=0, no=- 1, theta=0.7) Returns the gravitational acceleration field at position x by performing the Barnes-Hut treewalk using the provided octree instance Arguments: pos - (3,) array containing position of interest tree - octree instance storing the tree structure Keyword arguments: softening - softening radius of the particle at which the force is being evaluated - we use the greater of the target and source softenings when evaluating the softened potential no - index of the top-level node whose field is being summed - defaults to the global top-level node, can use a subnode in principle for e.g. parallelization theta - cell opening angle used to control force accuracy; smaller is slower (runtime ~ theta^-3) but more accurate. (default 0.7 gives ~1% accuracy)
-
pytreegrav.treewalk.
AccelWalk_quad
(pos, tree, softening=0, no=- 1, theta=0.7) Returns the gravitational acceleration field at position x by performing the Barnes-Hut treewalk using the provided octree instance. Uses the quadrupole expansion. Arguments: pos - (3,) array containing position of interest tree - octree instance storing the tree structure Keyword arguments: softening - softening radius of the particle at which the force is being evaluated - we use the greater of the target and source softenings when evaluating the softened potential no - index of the top-level node whose field is being summed - defaults to the global top-level node, can use a subnode in principle for e.g. parallelization theta - cell opening angle used to control force accuracy; smaller is slower (runtime ~ theta^-3) but more accurate. (default 0.7 gives ~1% accuracy)
-
pytreegrav.treewalk.
ColumnDensityWalk
(pos, rays, tree, no=- 1) Returns the integrated column density to infinity from pos, in the directions given by the rays argument
Arguments: pos - (3,) array containing position of interest rays - (N_rays, 3) array of unit vectors tree - octree object storing the tree structure
Returns: columns - (N_rays,) array of column densities along directions given by rays
Keyword arguments: no - index of the top-level node whose field is being summed - defaults to the global top-level node, can use a subnode in principle for e.g. parallelization
-
pytreegrav.treewalk.
ColumnDensity_tree
(pos_target, rays, tree) Returns the gravitational potential at the specified points, given a tree containing the mass distribution Arguments: pos_target – shape (N,3) array of positions at which to evaluate the potential tree – Octree instance containing the positions, masses, and softenings of the source particles
-
pytreegrav.treewalk.
ColumnDensity_tree_parallel
(pos_target, rays, tree) Returns the gravitational potential at the specified points, given a tree containing the mass distribution Arguments: pos_target – shape (N,3) array of positions at which to evaluate the potential tree – Octree instance containing the positions, masses, and softenings of the source particles
-
pytreegrav.treewalk.
DensityCorrFunc_tree
(pos, tree, rbins, max_bin_size_ratio=100, theta=0.7, boxsize=0, weighted_binning=False) Returns the average mass in radial bins surrounding a point
Arguments: pos – shape (N,3) array of particle positions tree – Octree instance containing the positions, masses, and softenings of the source particles
Optional arguments: rbins – 1D array of radial bin edges - if None will use heuristics to determine sensible bins max_bin_size_ratio – controls the accuracy of the binning - tree nodes are subdivided until their side length is at most this factor * the radial bin width
Returns: mbins – arrays containing total mass in each bin
-
pytreegrav.treewalk.
DensityCorrFunc_tree_parallel
(pos, tree, rbins, max_bin_size_ratio=100, theta=0.7, boxsize=0, weighted_binning=False) Returns the average mass in radial bins surrounding a point
Arguments: pos – shape (N,3) array of particle positions tree – Octree instance containing the positions, masses, and softenings of the source particles
Optional arguments: rbins – 1D array of radial bin edges - if None will use heuristics to determine sensible bins max_bin_size_ratio – controls the accuracy of the binning - tree nodes are subdivided until their side length is at most this factor * the radial bin width
Returns: mbins – arrays containing total mass in each bin
-
pytreegrav.treewalk.
DensityCorrWalk
(pos, tree, rbins, max_bin_size_ratio=100, theta=0.7, no=- 1, boxsize=0, weighted_binning=False) Returns the gravitational potential at position x by performing the Barnes-Hut treewalk using the provided octree instance
Arguments: pos - (3,) array containing position of interest tree - octree object storing the tree structure
Keyword arguments: softening - softening radius of the particle at which the force is being evaluated - we use the greater of the target and source softenings when evaluating the softened potential no - index of the top-level node whose field is being summed - defaults to the global top-level node, can use a subnode in principle for e.g. parallelization theta - cell opening angle used to control force accuracy; smaller is slower (runtime ~ theta^-3) but more accurate. (default 0.7 gives ~1% accuracy)
-
pytreegrav.treewalk.
NearestImage
(x, boxsize)
-
pytreegrav.treewalk.
PotentialTarget_tree
(pos_target, softening_target, tree, G=1.0, theta=0.7, quadrupole=False) Returns the gravitational potential at the specified points, given a tree containing the mass distribution Arguments: pos_target – shape (N,3) array of positions at which to evaluate the potential softening_target – shape (N,) array of minimum softening lengths to be used in all potential computations tree – Octree instance containing the positions, masses, and softenings of the source particles Optional arguments: G – gravitational constant (default 1.0) theta – accuracy parameter, smaller is more accurate, larger is faster (default 0.7) Returns: shape (N,) array of potential values at each point in pos
-
pytreegrav.treewalk.
PotentialTarget_tree_parallel
(pos_target, softening_target, tree, G=1.0, theta=0.7, quadrupole=False) Returns the gravitational potential at the specified points, given a tree containing the mass distribution Arguments: pos_target – shape (N,3) array of positions at which to evaluate the potential softening_target – shape (N,) array of minimum softening lengths to be used in all potential computations tree – Octree instance containing the positions, masses, and softenings of the source particles Optional arguments: G – gravitational constant (default 1.0) theta – accuracy parameter, smaller is more accurate, larger is faster (default 0.7) Returns: shape (N,) array of potential values at each point in pos
-
pytreegrav.treewalk.
PotentialWalk
(pos, tree, softening=0, no=- 1, theta=0.7) Returns the gravitational potential at position x by performing the Barnes-Hut treewalk using the provided octree instance Arguments: pos - (3,) array containing position of interest tree - octree object storing the tree structure Keyword arguments: softening - softening radius of the particle at which the force is being evaluated - we use the greater of the target and source softenings when evaluating the softened potential no - index of the top-level node whose field is being summed - defaults to the global top-level node, can use a subnode in principle for e.g. parallelization theta - cell opening angle used to control force accuracy; smaller is slower (runtime ~ theta^-3) but more accurate. (default 0.7 gives ~1% accuracy)
-
pytreegrav.treewalk.
PotentialWalk_quad
(pos, tree, softening=0, no=- 1, theta=0.7) Returns the gravitational potential at position x by performing the Barnes-Hut treewalk using the provided octree instance. Uses the quadrupole expansion. Arguments: pos - (3,) array containing position of interest tree - octree object storing the tree structure Keyword arguments: softening - softening radius of the particle at which the force is being evaluated - we use the greater of the target and source softenings when evaluating the softened potential no - index of the top-level node whose field is being summed - defaults to the global top-level node, can use a subnode in principle for e.g. parallelization theta - cell opening angle used to control force accuracy; smaller is slower (runtime ~ theta^-3) but more accurate. (default 0.7 gives ~1% accuracy)
-
pytreegrav.treewalk.
VelocityCorrFunc_tree
(pos, vel, weight, tree, rbins, max_bin_size_ratio=100, theta=0.7, boxsize=0, weighted_binning=False) Returns the average mass in radial bins surrounding a point
Arguments: pos – shape (N,3) array of particle positions tree – Octree instance containing the positions, masses, and softenings of the source particles
Optional arguments: rbins – 1D array of radial bin edges - if None will use heuristics to determine sensible bins max_bin_size_ratio – controls the accuracy of the binning - tree nodes are subdivided until their side length is at most this factor * the radial bin width (default 0.5)
Returns: mbins – arrays containing total mass in each bin
-
pytreegrav.treewalk.
VelocityCorrFunc_tree_parallel
(pos, vel, weight, tree, rbins, max_bin_size_ratio=100, theta=0.7, boxsize=0, weighted_binning=False) Returns the average mass in radial bins surrounding a point
Arguments: pos – shape (N,3) array of particle positions tree – Octree instance containing the positions, masses, and softenings of the source particles
Optional arguments: rbins – 1D array of radial bin edges - if None will use heuristics to determine sensible bins max_bin_size_ratio – controls the accuracy of the binning - tree nodes are subdivided until their side length is at most this factor * the radial bin width (default 0.5)
Returns: mbins – arrays containing total mass in each bin
-
pytreegrav.treewalk.
VelocityCorrWalk
(pos, vel, tree, rbins, max_bin_size_ratio=100, theta=0.7, no=- 1, boxsize=0, weighted_binning=False) Returns the gravitational potential at position x by performing the Barnes-Hut treewalk using the provided octree instance
Arguments: pos - (3,) array containing position of interest vel - (3,) array containing velocity of point of interest tree - octree object storing the tree structure
Keyword arguments: softening - softening radius of the particle at which the force is being evaluated - we use the greater of the target and source softenings when evaluating the softened potential no - index of the top-level node whose field is being summed - defaults to the global top-level node, can use a subnode in principle for e.g. parallelization theta - cell opening angle used to control force accuracy; smaller is slower (runtime ~ theta^-3) but more accurate. (default 0.7 gives ~1% accuracy)
-
pytreegrav.treewalk.
VelocityStructFunc_tree
(pos, vel, weight, tree, rbins, max_bin_size_ratio=100, theta=0.7, boxsize=0, weighted_binning=False) Returns the average mass in radial bins surrounding a point
Arguments: pos – shape (N,3) array of particle positions tree – Octree instance containing the positions, masses, and softenings of the source particles
Optional arguments: rbins – 1D array of radial bin edges - if None will use heuristics to determine sensible bins max_bin_size_ratio – controls the accuracy of the binning - tree nodes are subdivided until their side length is at most this factor * the radial bin width (default 0.5)
Returns: mbins – arrays containing total mass in each bin
-
pytreegrav.treewalk.
VelocityStructFunc_tree_parallel
(pos, vel, weight, tree, rbins, max_bin_size_ratio=100, theta=0.7, boxsize=0, weighted_binning=False) Returns the average mass in radial bins surrounding a point
Arguments: pos – shape (N,3) array of particle positions tree – Octree instance containing the positions, masses, and softenings of the source particles
Optional arguments: rbins – 1D array of radial bin edges - if None will use heuristics to determine sensible bins max_bin_size_ratio – controls the accuracy of the binning - tree nodes are subdivided until their side length is at most this factor * the radial bin width (default 0.5)
Returns: mbins – arrays containing total mass in each bin
-
pytreegrav.treewalk.
VelocityStructWalk
(pos, vel, tree, rbins, max_bin_size_ratio=100, theta=0.7, no=- 1, boxsize=0, weighted_binning=False) Returns the gravitational potential at position x by performing the Barnes-Hut treewalk using the provided octree instance
Arguments: pos - (3,) array containing position of interest vel - (3,) array containing velocity of point of interest tree - octree object storing the tree structure
Keyword arguments: softening - softening radius of the particle at which the force is being evaluated - we use the greater of the target and source softenings when evaluating the softened potential no - index of the top-level node whose field is being summed - defaults to the global top-level node, can use a subnode in principle for e.g. parallelization theta - cell opening angle used to control force accuracy; smaller is slower (runtime ~ theta^-3) but more accurate. (default 0.7 gives ~1% accuracy)
-
pytreegrav.treewalk.
do_weighted_binning
(tree, no, rbins, mbin, r, r_idx, quantity)