Module nnps: Nearest Neighbor Particle Search

class pysph.base.nnps.BoxSortNNPS

Bases: pysph.base.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

cell_to_index: ‘map[long,int]’

class pysph.base.nnps.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

centroid: ‘cPoint’

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

gindices: list

is_boundary

is_boundary: ‘bool’

lindices

lindices: list

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

Set the global and local indices for the cell

size

size: ‘int’

class pysph.base.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.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

cells: dict

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.nnps.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)

Bases: object

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

cell_size

cell_size: ‘double’

compute_cell_size_for_binning(self)
dim

dim: ‘int’

is_periodic

is_periodic: ‘bool’

narrays

narrays: ‘int’

pa_wrappers

pa_wrappers: list

periodic_in_x

periodic_in_x: ‘bool’

periodic_in_y

periodic_in_y: ‘bool’

periodic_in_z

periodic_in_z: ‘bool’

radius_scale

radius_scale: ‘double’

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)
update(self, *args, **kwargs)

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.

xmax

xmax: ‘double’

xmin

xmin: ‘double’

xtranslate

xtranslate: ‘double’

ymax

ymax: ‘double’

ymin

ymin: ‘double’

ytranslate

ytranslate: ‘double’

zmax

zmax: ‘double’

zmin

zmin: ‘double’

ztranslate

ztranslate: ‘double’

class pysph.base.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.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

fixed_h: ‘bool’

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.
get_spatially_ordered_indices(self, int pa_index, LongArray indices)
heads

heads: list

ncells_per_dim

ncells_per_dim: pyzoltan.core.carray.IntArray

ncells_tot

ncells_tot: ‘int’

nexts

nexts: list

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.NNPS(int dim, list particles, double radius_scale=2.0, int ghost_layers=1, domain=None, bool cache=False, bool sort_gids=False)

Bases: object

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.

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)
cell_size

cell_size: ‘double’

current_cache

current_cache: pysph.base.nnps.NeighborCache

dim

dim: ‘int’

domain

domain: pysph.base.nnps.DomainManager

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)
get_spatially_ordered_indices(self, int pa_index, LongArray indices)
is_periodic

is_periodic: ‘bool’

n_cells

n_cells: ‘int’

narrays

narrays: ‘int’

pa_wrappers

pa_wrappers: list

particles

particles: list

radius_scale

radius_scale: ‘double’

set_context(self, int src_index, int dst_index)

Setup the context before asing 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.) –
set_in_parallel(self, bool in_parallel)
sort_gids

sort_gids: ‘bool’

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, *args, **kwargs)
use_cache

use_cache: ‘bool’

xmax

xmax: pyzoltan.core.carray.DoubleArray

xmin

xmin: pyzoltan.core.carray.DoubleArray

class pysph.base.nnps.NNPSParticleArrayWrapper(ParticleArray pa)

Bases: object

gid

gid: pyzoltan.core.carray.UIntArray

h

h: pyzoltan.core.carray.DoubleArray

pa

pa: pysph.base.particle_array.ParticleArray

remove_tagged_particles(self, int tag)
tag

tag: pyzoltan.core.carray.IntArray

x

x: pyzoltan.core.carray.DoubleArray

y

y: pyzoltan.core.carray.DoubleArray

z

z: pyzoltan.core.carray.DoubleArray

class pysph.base.nnps.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.arange_uint(int start, int stop=-1) → UIntArray

Utility function to return a numpy.arange for a UIntArray

pysph.base.nnps.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.get_number_of_threads() → int
pysph.base.nnps.py_flatten(IntPoint cid, IntArray ncells_per_dim, int dim)

Python wrapper

pysph.base.nnps.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.py_unflatten(long cell_index, IntArray ncells_per_dim, int dim)

Python wrapper

pysph.base.nnps.set_number_of_threads(int n)