MFEM  v4.2.0
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
Public Types | Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
mfem::FindPointsGSLIB Class Reference

FindPointsGSLIB can robustly evaluate a GridFunction on an arbitrary collection of points. There are three key functions in FindPointsGSLIB: More...

#include <gslib.hpp>

Collaboration diagram for mfem::FindPointsGSLIB:
[legend]

Public Types

enum  AvgType { NONE, ARITHMETIC, HARMONIC }
 

Public Member Functions

 FindPointsGSLIB ()
 
 FindPointsGSLIB (MPI_Comm _comm)
 
 ~FindPointsGSLIB ()
 
void Setup (Mesh &m, const double bb_t=0.1, const double newt_tol=1.0e-12, const int npt_max=256)
 
void FindPoints (const Vector &point_pos)
 
void FindPoints (Mesh &m, const Vector &point_pos, const double bb_t=0.1, const double newt_tol=1.0e-12, const int npt_max=256)
 Setup FindPoints and search positions. More...
 
void Interpolate (const GridFunction &field_in, Vector &field_out)
 
void Interpolate (const Vector &point_pos, const GridFunction &field_in, Vector &field_out)
 
void Interpolate (Mesh &m, const Vector &point_pos, const GridFunction &field_in, Vector &field_out)
 
void SetL2AvgType (AvgType avgtype_)
 
void SetDefaultInterpolationValue (double interp_value_)
 
void FreeData ()
 
const Array< unsigned int > & GetCode () const
 
const Array< unsigned int > & GetElem () const
 Return element number for each point found by FindPoints. More...
 
const Array< unsigned int > & GetProc () const
 Return MPI rank on which each point was found by FindPoints. More...
 
const VectorGetReferencePosition () const
 Return reference coordinates for each point found by FindPoints. More...
 
const VectorGetDist () const
 
const Array< unsigned int > & GetGSLIBElem () const
 
const VectorGetGSLIBReferencePosition () const
 

Protected Member Functions

void GetNodeValues (const GridFunction &gf_in, Vector &node_vals)
 Get GridFunction from MFEM format to GSLIB format. More...
 
void GetQuadHexNodalCoordinates ()
 
void GetSimplexNodalCoordinates ()
 
void InterpolateH1 (const GridFunction &field_in, Vector &field_out)
 Use GSLIB for communication and interpolation. More...
 
void InterpolateGeneral (const GridFunction &field_in, Vector &field_out)
 
void MapRefPosAndElemIndices ()
 

Protected Attributes

Meshmesh
 
Meshmeshsplit
 
IntegrationRuleir_simplex
 
struct findpts_data_2 * fdata2D
 
struct findpts_data_3 * fdata3D
 
struct crystal * cr
 
struct comm * gsl_comm
 
int dim
 
int points_cnt
 
Array< unsigned int > gsl_code
 
Array< unsigned int > gsl_proc
 
Array< unsigned int > gsl_elem
 
Array< unsigned int > gsl_mfem_elem
 
Vector gsl_mesh
 
Vector gsl_ref
 
Vector gsl_dist
 
Vector gsl_mfem_ref
 
bool setupflag
 
double default_interp_value
 
AvgType avgtype
 

Detailed Description

FindPointsGSLIB can robustly evaluate a GridFunction on an arbitrary collection of points. There are three key functions in FindPointsGSLIB:

  1. Setup - constructs the internal data structures of gslib.
  2. FindPoints - for any given arbitrary set of points in physical space, gslib finds the element number, MPI rank, and the reference space coordinates inside the element that each point is located in. gslib also returns a code that indicates whether the point was found inside an element, on element border, or not found in the domain.
  3. Interpolate - Interpolates any grid function at the points found using 2.

FindPointsGSLIB provides interface to use these functions individually or using a single call.

Definition at line 45 of file gslib.hpp.

Member Enumeration Documentation

Enumerator
NONE 
ARITHMETIC 
HARMONIC 

Definition at line 48 of file gslib.hpp.

Constructor & Destructor Documentation

mfem::FindPointsGSLIB::FindPointsGSLIB ( )

Definition at line 31 of file gslib.cpp.

mfem::FindPointsGSLIB::FindPointsGSLIB ( MPI_Comm  _comm)

Definition at line 59 of file gslib.cpp.

mfem::FindPointsGSLIB::~FindPointsGSLIB ( )

Definition at line 50 of file gslib.cpp.

Member Function Documentation

void mfem::FindPointsGSLIB::FindPoints ( const Vector point_pos)

Searches positions given in physical space by point_pos. These positions must by ordered by nodes: (XXX...,YYY...,ZZZ). This function populates the following member variables: gsl_code Return codes for each point: inside element (0), element boundary (1), not found (2). gsl_proc MPI proc ids where the points were found. gsl_elem Element ids where the points were found. Defaults to 0 for points that were not found. gsl_mfem_elem Element ids corresponding to MFEM-mesh where the points were found. gsl_mfem_elem != gsl_elem for simplices Defaults to 0 for points that were not found. gsl_ref Reference coordinates of the found point. Ordered by vdim (XYZ,XYZ,XYZ...). Defaults to -1 for points that were not found. Note: the gslib reference frame is [-1,1]. gsl_mfem_ref Reference coordinates gsl_ref mapped to [0,1]. Defaults to 0 for points that were not found. gsl_dist Distance between the sought and the found point in physical space.

Definition at line 125 of file gslib.cpp.

void mfem::FindPointsGSLIB::FindPoints ( Mesh m,
const Vector point_pos,
const double  bb_t = 0.1,
const double  newt_tol = 1.0e-12,
const int  npt_max = 256 
)

Setup FindPoints and search positions.

Definition at line 183 of file gslib.cpp.

void mfem::FindPointsGSLIB::FreeData ( )

Cleans up memory allocated internally by gslib. Note that in parallel, this must be called before MPI_Finalize(), as it calls MPI_Comm_free() for internal gslib communicators.

Definition at line 208 of file gslib.cpp.

const Array<unsigned int>& mfem::FindPointsGSLIB::GetCode ( ) const
inline

Return code for each point searched by FindPoints: inside element (0), on element boundary (1), or not found (2).

Definition at line 162 of file gslib.hpp.

const Vector& mfem::FindPointsGSLIB::GetDist ( ) const
inline

Return distance Distance between the sought and the found point in physical space, for each point found by FindPoints.

Definition at line 171 of file gslib.hpp.

const Array<unsigned int>& mfem::FindPointsGSLIB::GetElem ( ) const
inline

Return element number for each point found by FindPoints.

Definition at line 164 of file gslib.hpp.

const Array<unsigned int>& mfem::FindPointsGSLIB::GetGSLIBElem ( ) const
inline

Return element number for each point found by FindPoints corresponding to GSLIB mesh. gsl_mfem_elem != gsl_elem for mesh with simplices.

Definition at line 175 of file gslib.hpp.

const Vector& mfem::FindPointsGSLIB::GetGSLIBReferencePosition ( ) const
inline

Return reference coordinates in [-1,1] (internal range in GSLIB) for each point found by FindPoints.

Definition at line 178 of file gslib.hpp.

void mfem::FindPointsGSLIB::GetNodeValues ( const GridFunction gf_in,
Vector node_vals 
)
protected

Get GridFunction from MFEM format to GSLIB format.

Definition at line 229 of file gslib.cpp.

const Array<unsigned int>& mfem::FindPointsGSLIB::GetProc ( ) const
inline

Return MPI rank on which each point was found by FindPoints.

Definition at line 166 of file gslib.hpp.

void mfem::FindPointsGSLIB::GetQuadHexNodalCoordinates ( )
protected

Get nodal coordinates from mesh to the format expected by GSLIB for quads and hexes

Definition at line 286 of file gslib.cpp.

const Vector& mfem::FindPointsGSLIB::GetReferencePosition ( ) const
inline

Return reference coordinates for each point found by FindPoints.

Definition at line 168 of file gslib.hpp.

void mfem::FindPointsGSLIB::GetSimplexNodalCoordinates ( )
protected

Convert simplices to quad/hexes and then get nodal coordinates for each split element into format expected by GSLIB

Definition at line 321 of file gslib.cpp.

void mfem::FindPointsGSLIB::Interpolate ( const GridFunction field_in,
Vector field_out 
)

Interpolation of field values at prescribed reference space positions.

Parameters
[in]field_inFunction values that will be interpolated on the reference positions. Note: it is assumed that field_in is in H1 and in the same space as the mesh that was given to Setup().
[out]field_outInterpolated values. For points that are not found the value is set to default_interp_value.

Definition at line 575 of file gslib.cpp.

void mfem::FindPointsGSLIB::Interpolate ( const Vector point_pos,
const GridFunction field_in,
Vector field_out 
)

Search positions and interpolate

Definition at line 194 of file gslib.cpp.

void mfem::FindPointsGSLIB::Interpolate ( Mesh m,
const Vector point_pos,
const GridFunction field_in,
Vector field_out 
)

Setup FindPoints, search positions and interpolate

Definition at line 201 of file gslib.cpp.

void mfem::FindPointsGSLIB::InterpolateGeneral ( const GridFunction field_in,
Vector field_out 
)
protected

Uses GSLIB Crystal Router for communication followed by MFEM's interpolation functions

Definition at line 692 of file gslib.cpp.

void mfem::FindPointsGSLIB::InterpolateH1 ( const GridFunction field_in,
Vector field_out 
)
protected

Use GSLIB for communication and interpolation.

Definition at line 650 of file gslib.cpp.

void mfem::FindPointsGSLIB::MapRefPosAndElemIndices ( )
protected

Map {r,s,t} coordinates from [-1,1] to [0,1] for MFEM. For simplices mesh find the original element number (that was split into micro quads/hexes by GetSimplexNodalCoordinates())

Definition at line 485 of file gslib.cpp.

void mfem::FindPointsGSLIB::SetDefaultInterpolationValue ( double  interp_value_)
inline

Set the default interpolation value for points that are not found in the mesh.

Definition at line 150 of file gslib.hpp.

void mfem::FindPointsGSLIB::SetL2AvgType ( AvgType  avgtype_)
inline

Average type to be used for L2 functions in-case a point is located at an element boundary where the function might be multi-valued.

Definition at line 146 of file gslib.hpp.

void mfem::FindPointsGSLIB::Setup ( Mesh m,
const double  bb_t = 0.1,
const double  newt_tol = 1.0e-12,
const int  npt_max = 256 
)

Initializes the internal mesh in gslib, by sending the positions of the Gauss-Lobatto nodes of the input Mesh object m. Note: not tested with periodic (DG meshes). Note: the input mesh m must have Nodes set.

Parameters
[in]mInput mesh.
[in]bb_tRelative size of bounding box around each element.
[in]newt_tolNewton tolerance for the gslib search methods.
[in]npt_maxNumber of points for simultaneous iteration. This alters performance and memory footprint.

Definition at line 71 of file gslib.cpp.

Member Data Documentation

AvgType mfem::FindPointsGSLIB::avgtype
protected

Definition at line 62 of file gslib.hpp.

struct crystal* mfem::FindPointsGSLIB::cr
protected

Definition at line 55 of file gslib.hpp.

double mfem::FindPointsGSLIB::default_interp_value
protected

Definition at line 61 of file gslib.hpp.

int mfem::FindPointsGSLIB::dim
protected

Definition at line 57 of file gslib.hpp.

struct findpts_data_2* mfem::FindPointsGSLIB::fdata2D
protected

Definition at line 53 of file gslib.hpp.

struct findpts_data_3* mfem::FindPointsGSLIB::fdata3D
protected

Definition at line 54 of file gslib.hpp.

Array<unsigned int> mfem::FindPointsGSLIB::gsl_code
protected

Definition at line 58 of file gslib.hpp.

struct comm* mfem::FindPointsGSLIB::gsl_comm
protected

Definition at line 56 of file gslib.hpp.

Vector mfem::FindPointsGSLIB::gsl_dist
protected

Definition at line 59 of file gslib.hpp.

Array<unsigned int> mfem::FindPointsGSLIB::gsl_elem
protected

Definition at line 58 of file gslib.hpp.

Vector mfem::FindPointsGSLIB::gsl_mesh
protected

Definition at line 59 of file gslib.hpp.

Array<unsigned int> mfem::FindPointsGSLIB::gsl_mfem_elem
protected

Definition at line 58 of file gslib.hpp.

Vector mfem::FindPointsGSLIB::gsl_mfem_ref
protected

Definition at line 59 of file gslib.hpp.

Array<unsigned int> mfem::FindPointsGSLIB::gsl_proc
protected

Definition at line 58 of file gslib.hpp.

Vector mfem::FindPointsGSLIB::gsl_ref
protected

Definition at line 59 of file gslib.hpp.

IntegrationRule* mfem::FindPointsGSLIB::ir_simplex
protected

Definition at line 52 of file gslib.hpp.

Mesh* mfem::FindPointsGSLIB::mesh
protected

Definition at line 51 of file gslib.hpp.

Mesh * mfem::FindPointsGSLIB::meshsplit
protected

Definition at line 51 of file gslib.hpp.

int mfem::FindPointsGSLIB::points_cnt
protected

Definition at line 57 of file gslib.hpp.

bool mfem::FindPointsGSLIB::setupflag
protected

Definition at line 60 of file gslib.hpp.


The documentation for this class was generated from the following files: