20 #ifndef LIBMESH_MESHFREE_INTERPOLATION_H 21 #define LIBMESH_MESHFREE_INTERPOLATION_H 24 #include "libmesh/libmesh_config.h" 25 #include "libmesh/libmesh_common.h" 26 #include "libmesh/point.h" 27 #include "libmesh/parallel_object.h" 28 #ifdef LIBMESH_HAVE_NANOFLANN 29 # include "libmesh/ignore_warnings.h" 30 # include "libmesh/nanoflann.hpp" 31 # include "libmesh/restore_warnings.h" 91 friend std::ostream &
operator << (std::ostream & os,
104 {
return cast_int<unsigned int>(
_names.size()); }
111 {
_names = std::move(names); }
134 virtual void add_field_data (
const std::vector<std::string> & field_names,
135 const std::vector<Point> & pts,
136 const std::vector<Number> & vals);
151 const std::vector<Point> & tgt_pts,
152 std::vector<Number> & tgt_vals)
const = 0;
178 template <
unsigned int KDDim>
183 #ifdef LIBMESH_HAVE_NANOFLANN 190 template <
unsigned int PLDim>
194 const std::vector<Point> &
_pts;
217 libmesh_assert_equal_to (size, PLDim);
218 libmesh_assert_less (idx_p2,
_pts.size());
226 const coord_t d0=p1[0] - p2(0);
227 const coord_t d1=p1[1] - p2(1);
228 const coord_t d2=p1[2] - p2(2);
230 return d0*d0 + d1*d1 + d2*d2;
235 const coord_t d0=p1[0] - p2(0);
236 const coord_t d1=p1[1] - p2(1);
238 return d0*d0 + d1*d1;
243 const coord_t d0=p1[0] - p2(0);
249 libmesh_error_msg(
"ERROR: unknown size " << size);
262 libmesh_assert_less (
dim, PLDim);
263 libmesh_assert_less (
idx,
_pts.size());
264 libmesh_assert_less (
dim, 3);
268 if (
dim==0)
return p(0);
269 if (
dim==1)
return p(1);
283 template <
class BBOX>
296 typedef nanoflann::KDTreeSingleIndexAdaptor<nanoflann::L2_Simple_Adaptor<Real, PointListAdaptor<KDDim>>,
301 #endif // LIBMESH_HAVE_NANOFLANN 313 const std::vector<size_t> & src_indices,
314 const std::vector<Real> & src_dist_sqr,
315 std::vector<Number>::iterator & out_it)
const;
332 const unsigned int n_interp_pts = 8,
333 const Real power = 2) :
335 #
if LIBMESH_HAVE_NANOFLANN
346 virtual void clear()
override;
353 const std::vector<Point> & tgt_pts,
354 std::vector<Number> & tgt_vals)
const override;
360 #endif // #define LIBMESH_MESHFREE_INTERPOLATION_H nanoflann::KDTreeSingleIndexAdaptor< nanoflann::L2_Simple_Adaptor< Real, PointListAdaptor< KDDim > >, PointListAdaptor< KDDim >, KDDim, std::size_t > kd_tree_t
unsigned int n_field_variables() const
The number of field variables.
Inverse distance interpolation.
const std::vector< std::string > & field_variables() const
std::vector< Number > _vals
Temporary work array.
size_t kdtree_get_point_count() const
Must return the number of data points.
The libMesh namespace provides an interface to certain functionality in the library.
coord_t kdtree_distance(const coord_t *p1, const size_t idx_p2, size_t size) const
virtual void interpolate(const Point &pt, const std::vector< size_t > &src_indices, const std::vector< Real > &src_dist_sqr, std::vector< Number >::iterator &out_it) const
Performs inverse distance interpolation at the input point from the specified points.
coord_t kdtree_get_pt(const size_t idx, int dim) const
std::vector< Point > & get_source_points()
const std::vector< Point > & _pts
virtual void interpolate_field_data(const std::vector< std::string > &field_names, const std::vector< Point > &tgt_pts, std::vector< Number > &tgt_vals) const =0
Interpolate source data at target points.
virtual void clear()
Clears all internal data structures and restores to a pristine state.
ParallelizationStrategy _parallelization_strategy
std::vector< Number > _src_vals
virtual void add_field_data(const std::vector< std::string > &field_names, const std::vector< Point > &pts, const std::vector< Number > &vals)
Sets source data at specified points.
PointListAdaptor(const std::vector< Point > &pts)
std::vector< Point > _src_pts
PointListAdaptor< KDDim > _point_list_adaptor
void print_info(std::ostream &os=libMesh::out) const
Prints information about this object, by default to libMesh::out.
Base class to support various mesh-free interpolation methods.
virtual void construct_kd_tree()
Build & initialize the KD tree, if needed.
virtual void clear() override
Clears all internal data structures and restores to a pristine state.
friend std::ostream & operator<<(std::ostream &os, const MeshfreeInterpolation &mfi)
Same as above, but allows you to also use stream syntax.
std::vector< Number > & get_source_vals()
This class adapts list of libMesh Point types for use in a nanoflann KD-Tree.
virtual void gather_remote_data()
Gathers source points and values that have been added on other processors.
An object whose state is distributed along a set of processors.
virtual void prepare_for_use()
Prepares data structures for use.
virtual void interpolate_field_data(const std::vector< std::string > &field_names, const std::vector< Point > &tgt_pts, std::vector< Number > &tgt_vals) const override
Interpolate source data at target points.
InverseDistanceInterpolation(const libMesh::Parallel::Communicator &comm_in, const unsigned int n_interp_pts=8, const Real power=2)
Constructor.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void set_field_variables(std::vector< std::string > names)
Defines the field variable(s) we are responsible for, and importantly their assumed ordering...
Real coord_t
libMesh Point coordinate type
std::vector< std::string > _names
std::unique_ptr< kd_tree_t > _kd_tree
MeshfreeInterpolation(const libMesh::Parallel::Communicator &comm_in)
Constructor.
ParallelizationStrategy
"ParallelizationStrategy" to employ.
A Point defines a location in LIBMESH_DIM dimensional Real space.
const unsigned int _n_interp_pts
bool kdtree_get_bbox(BBOX &) const
Optional bounding-box computation.