Module nnps: Nearest Neighbor Particle Search

class pysph.base.nnps_base.CPUDomainManager(double xmin=-1000, double xmax=1000, double ymin=0, double ymax=0, double zmin=0, double zmax=0, periodic_in_x=False, periodic_in_y=False, periodic_in_z=False, double n_layers=2.0, backend=None, props=None, mirror_in_x=False, mirror_in_y=False, mirror_in_z=False)

Bases: pysph.base.nnps_base.DomainManagerBase

This class determines the limits of the solution domain.

We expect all simulations to have well defined domain limits beyond which we are either not interested or the solution is invalid to begin with. Thus, if a particle leaves the domain, the solution should be considered invalid (at least locally).

The initial domain limits could be given explicitly or asked to be computed from the particle arrays. The domain could be periodic.

Constructor

The n_layers argument specifies the number of ghost layers as multiples of hmax*radius_scale.

props: list/dict: properties to copy.
Provide a list or dict with the keys as particle array names. Only the specified properties are copied. If not specified, all props are copied.
dtype

object

Type:dtype
dtype_max

‘double’

Type:dtype_max
ghosts

list

Type:ghosts
update(self)

General method that is called before NNPS can bin particles.

This method is responsible for the computation of cell sizes and creation of any ghost particles for periodic or wall boundary conditions.

use_double

‘bool’

Type:use_double
class pysph.base.nnps_base.Cell(IntPoint cid, double cell_size, int narrays, int layers=2)

Bases: object

Basic indexing structure for the box-sort NNPS.

For a spatial indexing based on the box-sort algorithm, this class defines the spatial data structure used to hold particle indices (local and global) that are within this cell.

Constructor

Parameters:
  • cid (IntPoint) – Spatial index (unflattened) for the cell
  • cell_size (double) – Spatial extent of the cell in each dimension
  • narrays (int) – Number of arrays being binned
  • layers (int) – Factor to compute the bounding box
centroid

‘cPoint’

Type:centroid
get_bounding_box(self, Point boxmin, Point boxmax, int layers=1, cell_size=None)

Compute the bounding box for the cell.

Parameters:
  • boxmin (Point (output)) – The bounding box min coordinates are stored here
  • boxmax (Point (output)) – The bounding box max coordinates are stored here
  • layers (int (input) default (1)) – Number of offset layers to define the bounding box
  • cell_size (double (input) default (None)) – Optional cell size to use to compute the bounding box. If not provided, the cell’s size will be used.
get_centroid(self, Point pnt)

Utility function to get the centroid of the cell.

Parameters:
  • pnt (Point (input/output)) – The centroid is cmoputed and stored in this object.
  • centroid is defined as the origin plus half the cell size (The) –
  • each dimension. (in) –
gindices

list

Type:gindices
is_boundary

‘bool’

Type:is_boundary
lindices

list

Type:lindices
set_indices(self, int index, UIntArray lindices, UIntArray gindices)

Set the global and local indices for the cell

size

‘double’

Type:size
class pysph.base.nnps_base.DomainManager(double xmin=-1000, double xmax=1000, double ymin=0, double ymax=0, double zmin=0, double zmax=0, periodic_in_x=False, periodic_in_y=False, periodic_in_z=False, double n_layers=2.0, backend=None, props=None, mirror_in_x=False, mirror_in_y=False, mirror_in_z=False)

Bases: object

Constructor

Parameters:
  • xmax, ymin, ymax, zmin, zmax (xmin,) –
  • periodic_in_y, periodic_in_z (periodic_in_x,) –
  • mirror_in_y, mirror_in_z (mirror_in_x,) –
  • n_layers (double: number of ghost layers as a multiple of) – h_max*radius_scale
  • props (list/dict: properties to copy.) – Provide a list or dict with the keys as particle array names. Only the specified properties are copied. If not specified, all props are copied.
backend

object

Type:backend
compute_cell_size_for_binning(self)

Compute the cell size for the binning.

The cell size is chosen as the kernel radius scale times the maximum smoothing length in the local processor. For parallel runs, we would need to communicate the maximum ‘h’ on all processors to decide on the appropriate binning size.

manager

object

Type:manager
set_cell_size(self, cell_size)
set_in_parallel(self, in_parallel)
set_pa_wrappers(self, wrappers)
set_radius_scale(self, radius_scale)
update(self)
class pysph.base.nnps_base.DomainManagerBase(double xmin=-1000, double xmax=1000, double ymin=0, double ymax=0, double zmin=0, double zmax=0, periodic_in_x=False, periodic_in_y=False, periodic_in_z=False, double n_layers=2.0, props=None, mirror_in_x=False, mirror_in_y=False, mirror_in_z=False)

Bases: object

Constructor

The n_layers argument specifies the number of ghost layers as multiples of hmax*radius_scale.

props: list/dict: properties to copy.
Provide a list or dict with the keys as particle array names. Only the specified properties are copied. If not specified, all props are copied.
cell_size

‘double’

Type:cell_size
compute_cell_size_for_binning(self)
copy_props

list

Type:copy_props
dim

‘int’

Type:dim
hmin

‘double’

Type:hmin
in_parallel

‘bool’

Type:in_parallel
is_mirror

‘bool’

Type:is_mirror
is_periodic

‘bool’

Type:is_periodic
mirror_in_x

‘bool’

Type:mirror_in_x
mirror_in_y

‘bool’

Type:mirror_in_y
mirror_in_z

‘bool’

Type:mirror_in_z
n_layers

‘double’

Type:n_layers
narrays

‘int’

Type:narrays
pa_wrappers

list

Type:pa_wrappers
periodic_in_x

‘bool’

Type:periodic_in_x
periodic_in_y

‘bool’

Type:periodic_in_y
periodic_in_z

‘bool’

Type:periodic_in_z
props

object

Type:props
radius_scale

‘double’

Type:radius_scale
set_cell_size(self, cell_size)
set_in_parallel(self, bool in_parallel)
set_pa_wrappers(self, wrappers)
set_radius_scale(self, double radius_scale)
xmax

‘double’

Type:xmax
xmin

‘double’

Type:xmin
xtranslate

‘double’

Type:xtranslate
ymax

‘double’

Type:ymax
ymin

‘double’

Type:ymin
ytranslate

‘double’

Type:ytranslate
zmax

‘double’

Type:zmax
zmin

‘double’

Type:zmin
ztranslate

‘double’

Type:ztranslate
class pysph.base.nnps_base.NNPS(int dim, list particles, double radius_scale=2.0, int ghost_layers=1, domain=None, bool cache=False, bool sort_gids=False)

Bases: pysph.base.nnps_base.NNPSBase

Nearest neighbor query class using the box-sort algorithm.

NNPS bins all local particles using the box sort algorithm in Cells. The cells are stored in a dictionary ‘cells’ which is keyed on the spatial index (IntPoint) of the cell.

current_cache

pysph.base.nnps_base.NeighborCache

Type:current_cache
get_spatially_ordered_indices(self, int pa_index, LongArray indices)
set_in_parallel(self, bool in_parallel)
set_use_cache(self, bool use_cache)
sort_gids

‘bool’

Type:sort_gids
spatially_order_particles(self, int pa_index)

Spatially order particles such that nearby particles have indices nearer each other. This may improve pre-fetching on the CPU.

update(self)

Update the local data after particles have moved.

For parallel runs, we want the NNPS to be independent of the ParallelManager which is solely responsible for distributing particles across available processors. We assume therefore that after a parallel update, each processor has all the local particle information it needs and this operation is carried out locally.

For serial runs, this method should be called when the particles have moved.

update_domain(self)
xmax

cyarray.carray.DoubleArray

Type:xmax
xmin

cyarray.carray.DoubleArray

Type:xmin
class pysph.base.nnps_base.NNPSBase(int dim, list particles, double radius_scale=2.0, int ghost_layers=1, domain=None, bool cache=False, bool sort_gids=False)

Bases: object

Constructor for NNPS

Parameters:
  • dim (int) – Dimension (fixme: Not sure if this is really needed)
  • particles (list) – The list of particle arrays we are working on.
  • radius_scale (double, default (2)) – Optional kernel radius scale. Defaults to 2
  • domain (DomainManager, default (None)) – Optional limits for the domain
  • cache (bint) – Flag to set if we want to cache neighbor calls. This costs storage but speeds up neighbor calculations.
  • sort_gids (bint, default (False)) – Flag to sort neighbors using gids (if they are available). This is useful when comparing parallel results with those from a serial run.
brute_force_neighbors(self, int src_index, int dst_index, size_t d_idx, UIntArray nbrs)
cache

list

Type:cache
cell_size

‘double’

Type:cell_size
dim

‘int’

Type:dim
domain

pysph.base.nnps_base.DomainManager

Type:domain
get_nearest_particles(self, int src_index, int dst_index, size_t d_idx, UIntArray nbrs)
get_nearest_particles_no_cache(self, int src_index, int dst_index, size_t d_idx, UIntArray nbrs, bool prealloc)

Find nearest neighbors for particle id ‘d_idx’ without cache

Parameters:
  • src_index (int) – Index in the list of particle arrays to which the neighbors belong
  • dst_index (int) – Index in the list of particle arrays to which the query point belongs
  • d_idx (size_t) – Index of the query point in the destination particle array
  • nbrs (UIntArray) – Array to be populated by nearest neighbors of ‘d_idx’
hmin

‘double’

Type:hmin
is_periodic

‘bool’

Type:is_periodic
n_cells

‘int’

Type:n_cells
narrays

‘int’

Type:narrays
pa_wrappers

list

Type:pa_wrappers
particles

list

Type:particles
radius_scale

‘double’

Type:radius_scale
set_context(self, int src_index, int dst_index)

Setup the context before asking for neighbors. The dst_index represents the particles for whom the neighbors are to be determined from the particle array with index src_index.

Parameters:
  • src_index (int: the source index of the particle array.) –
  • dst_index (int: the destination index of the particle array.) –
class pysph.base.nnps_base.NNPSParticleArrayWrapper(ParticleArray pa)

Bases: object

gid

cyarray.carray.UIntArray

Type:gid
h

cyarray.carray.DoubleArray

Type:h
name

unicode

Type:name
np

‘int’

Type:np
pa

pysph.base.particle_array.ParticleArray

Type:pa
remove_tagged_particles(self, int tag)
tag

cyarray.carray.IntArray

Type:tag
x

cyarray.carray.DoubleArray

Type:x
y

cyarray.carray.DoubleArray

Type:y
z

cyarray.carray.DoubleArray

Type:z
class pysph.base.nnps_base.NeighborCache(NNPS nnps, int dst_index, int src_index)

Bases: object

find_all_neighbors(self)
get_neighbors(self, int src_index, size_t d_idx, UIntArray nbrs)
update(self)
pysph.base.nnps_base.arange_uint(int start, int stop=-1) → UIntArray

Utility function to return a numpy.arange for a UIntArray

pysph.base.nnps_base.get_centroid(double cell_size, IntPoint cid)

Get the centroid of the cell.

Parameters:
  • cell_size (double (input)) – Cell size used for binning
  • cid (IntPoint (input)) – Spatial index for a cell
Returns:

centroid

Return type:

Point

Notes

The centroid in any coordinate direction is defined to be the origin plus half the cell size in that direction

pysph.base.nnps_base.get_number_of_threads() → int
pysph.base.nnps_base.py_flatten(IntPoint cid, IntArray ncells_per_dim, int dim)

Python wrapper

pysph.base.nnps_base.py_get_valid_cell_index(IntPoint cid, IntArray ncells_per_dim, int dim, int n_cells)

Return the flattened cell index for a valid cell

pysph.base.nnps_base.py_unflatten(long cell_index, IntArray ncells_per_dim, int dim)

Python wrapper

pysph.base.nnps_base.set_number_of_threads(int n)
class pysph.base.linked_list_nnps.LinkedListNNPS(int dim, list particles, double radius_scale=2.0, int ghost_layers=1, domain=None, bool fixed_h=False, bool cache=False, bool sort_gids=False)

Bases: pysph.base.nnps_base.NNPS

Nearest neighbor query class using the linked list method.

Constructor for NNPS

Parameters:
  • dim (int) – Number of dimension.
  • particles (list) – The list of particle arrays we are working on
  • radius_scale (double, default (2)) – Optional kernel radius scale. Defaults to 2
  • ghost_layers (int) – Optional number of layers to share in parallel
  • domain (DomainManager, default (None)) – Optional limits for the domain
  • fixed_h (bint) – Optional flag to use constant cell sizes throughout.
  • cache (bint) – Flag to set if we want to cache neighbor calls. This costs storage but speeds up neighbor calculations.
  • sort_gids (bint, default (False)) – Flag to sort neighbors using gids (if they are available). This is useful when comparing parallel results with those from a serial run.
fixed_h

‘bool’

Type:fixed_h
get_spatially_ordered_indices(self, int pa_index, LongArray indices)
heads

list

Type:heads
ncells_per_dim

cyarray.carray.IntArray

Type:ncells_per_dim
ncells_tot

‘int’

Type:ncells_tot
nexts

list

Type:nexts
set_context(self, int src_index, int dst_index)

Setup the context before asking for neighbors. The dst_index represents the particles for whom the neighbors are to be determined from the particle array with index src_index.

Parameters:
  • src_index (int: the source index of the particle array.) –
  • dst_index (int: the destination index of the particle array.) –
class pysph.base.box_sort_nnps.BoxSortNNPS

Bases: pysph.base.linked_list_nnps.LinkedListNNPS

Nearest neighbor query class using the box sort method but which uses the LinkedList algorithm. This makes this very fast but perhaps not as safe as the DictBoxSortNNPS. All this class does is to use a std::map to obtain a linear cell index from the actual flattened cell index.

Constructor for NNPS

Parameters:
  • dim (int) – Number of dimension.
  • particles (list) – The list of particle arrays we are working on
  • radius_scale (double, default (2)) – Optional kernel radius scale. Defaults to 2
  • ghost_layers (int) – Optional number of layers to share in parallel
  • domain (DomainManager, default (None)) – Optional limits for the domain
  • fixed_h (bint) – Optional flag to use constant cell sizes throughout.
  • cache (bint) – Flag to set if we want to cache neighbor calls. This costs storage but speeds up neighbor calculations.
  • sort_gids (bint, default (False)) – Flag to sort neighbors using gids (if they are available). This is useful when comparing parallel results with those from a serial run.
cell_to_index

‘map[long,int]’

Type:cell_to_index
class pysph.base.box_sort_nnps.DictBoxSortNNPS(int dim, list particles, double radius_scale=2.0, int ghost_layers=1, domain=None, cache=False, sort_gids=False)

Bases: pysph.base.nnps_base.NNPS

Nearest neighbor query class using the box-sort algorithm using a dictionary.

NNPS bins all local particles using the box sort algorithm in Cells. The cells are stored in a dictionary ‘cells’ which is keyed on the spatial index (IntPoint) of the cell.

Constructor for NNPS

Parameters:
  • dim (int) – Number of dimensions.
  • particles (list) – The list of particle arrays we are working on.
  • radius_scale (double, default (2)) – Optional kernel radius scale. Defaults to 2
  • domain (DomainManager, default (None)) – Optional limits for the domain
  • cache (bint) – Flag to set if we want to cache neighbor calls. This costs storage but speeds up neighbor calculations.
  • sort_gids (bint, default (False)) – Flag to sort neighbors using gids (if they are available). This is useful when comparing parallel results with those from a serial run.
cells

dict

Type:cells
get_nearest_particles_no_cache(self, int src_index, int dst_index, size_t d_idx, UIntArray nbrs, bool prealloc)

Utility function to get near-neighbors for a particle.

Parameters:
  • src_index (int) – Index of the source particle array in the particles list
  • dst_index (int) – Index of the destination particle array in the particles list
  • d_idx (int (input)) – Destination particle for which neighbors are sought.
  • nbrs (UIntArray (output)) – Neighbors for the requested particle are stored here.
  • prealloc (bool) – Specifies if the neighbor array already has pre-allocated space for the neighbor list. In this case the neighbors are directly set in the given array without resetting or appending to the array. This improves performance when the neighbors are cached.
class pysph.base.spatial_hash_nnps.ExtendedSpatialHashNNPS(int dim, list particles, double radius_scale=2.0, int H=3, int ghost_layers=1, domain=None, bool fixed_h=False, bool cache=False, bool sort_gids=False, long long table_size=131072, bool approximate=False)

Bases: pysph.base.nnps_base.NNPS

Finds nearest neighbors using Extended Spatial Hashing algorithm

Sub-divides each cell into smaller ones. Useful when particles cluster in a cell.

For approximate Extended Spatial Hash, if the distance between a cell and the cell of the query particle is greater than search radius, the entire cell is ignored.

Ref. http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.105.6732&rep=rep1&type=pdf

set_context(self, int src_index, int dst_index)

Set context for nearest neighbor searches.

Parameters:
  • src_index (int) – Index in the list of particle arrays to which the neighbors belong
  • dst_index (int) – Index in the list of particle arrays to which the query point belongs
class pysph.base.spatial_hash_nnps.SpatialHashNNPS(int dim, list particles, double radius_scale=2.0, int ghost_layers=1, domain=None, bool fixed_h=False, bool cache=False, bool sort_gids=False, long long table_size=131072)

Bases: pysph.base.nnps_base.NNPS

Nearest neighbor particle search using Spatial Hashing algorithm

Uses a hashtable to store particles according to cell it belongs to.

Ref. http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.105.6732&rep=rep1&type=pdf

set_context(self, int src_index, int dst_index)

Set context for nearest neighbor searches.

Parameters:
  • src_index (int) – Index in the list of particle arrays to which the neighbors belong
  • dst_index (int) – Index in the list of particle arrays to which the query point belongs