DMPlexFindVertices#
Try to find DAG points based on their coordinates.
Synopsis#
#include "petscdmplex.h"
#include "petscfe.h"
PetscErrorCode DMPlexFindVertices(DM dm, Vec coordinates, PetscReal eps, IS *points)
Not Collective (provided DMGetCoordinatesLocalSetUp() has been already called)
Input Parameters#
dm - The
DMPLEXobjectcoordinates - The
Vecof coordinates of the sought pointseps - The tolerance or
PETSC_DEFAULT
Output Parameter#
points - The
ISof found DAG points or -1
Notes#
The length of Vec coordinates must be npoints * dim where dim is the spatial dimension returned by DMGetCoordinateDim() and npoints is the number of sought points.
The output IS is living on PETSC_COMM_SELF and its length is npoints.
Each rank does the search independently.
If this rank’s local DMPLEX portion contains the DAG point corresponding to the i-th tuple of coordinates, the i-th entry of the output IS is set to that DAG point, otherwise to -1.
The output IS must be destroyed by user.
The tolerance is interpreted as the maximum Euclidean (L2) distance of the sought point from the specified coordinates.
Complexity of this function is currently O(mn) with m number of vertices to find and n number of vertices in the local mesh. This could probably be improved if needed.
See Also#
Level#
intermediate
Location#
src/dm/impls/plex/plexgeometry.c
Index of all DMPlex routines
Table of Contents for all manual pages
Index of all manual pages