///////////////////////////////////////////////////////////////////// ///////////////////////////////////////////////////////////////////// More...
#include <sample_point_container.h>
Public Member Functions | |
RefineableBinArray (SamplePointContainerParameters *bin_array_parameters_pt) | |
Constructor. More... | |
RefineableBinArray (const RefineableBinArray &data)=delete | |
Broken copy constructor. More... | |
void | operator= (const RefineableBinArray &)=delete |
Broken assignment operator. More... | |
~RefineableBinArray () | |
Destructor. More... | |
RefineableBinArray * | root_bin_array_pt () const |
Root bin array. More... | |
RefineableBin * | bin_pt (const unsigned &i) const |
Pointer to i-th bin; can be null if bin is empty. More... | |
unsigned | nbin () const |
Number of bins (not taking recursion into account) More... | |
unsigned | total_number_of_sample_points_computed_recursively () const |
Compute total number of sample points recursively. More... | |
void | fill_bin_array (const Vector< SamplePoint * > &sample_point_pt) |
Fill the bin array with specified SamplePoints. More... | |
void | add_sample_point (SamplePoint *new_sample_point_pt, const Vector< double > &zeta) |
Add specified SamplePoint to RefineableBinArray. More... | |
void | locate_zeta (const Vector< double > &zeta, GeomObject *&sub_geom_object_pt, Vector< double > &s) |
Find sub-GeomObject (finite element) and the local coordinate s within it that contains point with global coordinate zeta. sub_geom_object_pt=0 if point can't be found. More... | |
void | get_bin_boundaries (const unsigned &bin_index, Vector< std::pair< double, double >> &min_and_max_coordinates) |
Boundaries of specified bin in each coordinate direction. *.first = min; *.second = max. More... | |
unsigned | depth () const |
Depth of the hierarchical bin_array. More... | |
unsigned | max_depth () const |
Max depth of the hierarchical bin_array; const version. More... | |
unsigned & | max_depth () |
Max depth of the hierarchical bin_array. More... | |
bool | bin_array_is_recursive () const |
Is the BinArray recursive? More... | |
unsigned | max_number_of_sample_point_per_bin () const |
Maximum number of sample points in bin (before its subdivided recursively) More... | |
void | output_bins (std::ofstream &outfile) |
Output bins. More... | |
void | output_bin_vertices (std::ofstream &outfile) |
Output bin vertices (allowing display of bins as zones). More... | |
void | output_neighbouring_bins (const unsigned &bin_index, const unsigned &radius, std::ofstream &outfile) |
Output neighbouring bins at given "radius" of the specified bin. More... | |
unsigned & | total_number_of_sample_points_visited_during_locate_zeta_from_top_level () |
Counter to keep track of how many sample points we've visited during top level call to locate_zeta. More... | |
unsigned & | first_sample_point_to_actually_lookup_during_locate_zeta () |
When searching through sample points recursively from the top level RefineableBinArray (in deterministic order!) only actually do the locate_zeta calls when when the counter exceeds this value. More... | |
unsigned & | last_sample_point_to_actually_lookup_during_locate_zeta () |
When searching through sample points recursively from the top level RefineableBinArray (in deterministic order!) only actually do the locate_zeta calls when when the counter is less than this value. More... | |
unsigned & | multiplier_for_max_sample_point_to_actually_lookup_during_locate_zeta () |
Every time we've completed a "spiral", visiting a finite number of sample points in a deterministic order, use this multiplier to increase the max. number of sample points to be visited. Using a multiplier rather than a constant increment increases the amount of (more and more unlikely to yield anything!) work done locally before doing another costly mpi round trip when we're already far from the point we're trying to find. More... | |
unsigned & | initial_last_sample_point_to_actually_lookup_during_locate_zeta () |
When searching through sample points recursively from the top level RefineableBinArray (in deterministic order!) only actually do the locate_zeta calls when when the counter exceeds this value. This is the initial value when starting the spiral based search. More... | |
Public Member Functions inherited from BinArray | |
BinArray (Mesh *mesh_pt, const Vector< std::pair< double, double >> &min_and_max_coordinates, const Vector< unsigned > &dimensions_of_bin_array, const bool &use_eulerian_coordinates_during_setup, const bool &ignore_halo_elements_during_locate_zeta_search, const unsigned &nsample_points_generated_per_element) | |
Constructor. More... | |
BinArray () | |
Broken default constructor; needed for broken copy constructors. Don't call. It will die. More... | |
BinArray (const BinArray &data)=delete | |
Broken copy constructor. More... | |
void | operator= (const BinArray &)=delete |
Broken assignment operator. More... | |
virtual | ~BinArray () |
Virtual destructor. More... | |
void | get_neighbouring_bins_helper (const unsigned &bin_index, const unsigned &radius, Vector< unsigned > &neighbouring_bin_index, const bool &use_old_version=true) |
Helper function for computing the bin indices of neighbouring bins at a given "radius" of the specified bin. Final, optional boolean (default: true) chooses to use the old version which appears to be faster than Louis' new one after all (in the few cases where this functionality is still used – not all if we have cgal!) More... | |
void | profile_get_neighbouring_bins_helper () |
Profiling function to compare performance of two different versions of the get_neighbouring_bins_helper(...) function. More... | |
unsigned | coords_to_bin_index (const Vector< double > &zeta) |
Get (linearly enumerated) bin index of bin that contains specified zeta. More... | |
void | coords_to_vectorial_bin_index (const Vector< double > &zeta, Vector< unsigned > &bin_index) |
Get "coordinates" of bin that contains specified zeta. More... | |
unsigned | max_bin_dimension () const |
Max. bin dimension (number of bins in coordinate directions) More... | |
unsigned | ndim_zeta () const |
Dimension of the zeta ( = dim of local coordinate of elements) More... | |
unsigned | dimension_of_bin_array (const unsigned &i) const |
Number of bins in coordinate direction i. More... | |
Vector< unsigned > | dimensions_of_bin_array () const |
Number of bins in coordinate directions. Const vector-based version. More... | |
unsigned | dimensions_of_bin_array (const unsigned &i) const |
Number of bins in specified coordinate direction. More... | |
Public Member Functions inherited from SamplePointContainer | |
SamplePointContainer (Mesh *mesh_pt, const Vector< std::pair< double, double >> &min_and_max_coordinates, const bool &use_eulerian_coordinates_during_setup, const bool &ignore_halo_elements_during_locate_zeta_search, const unsigned &nsample_points_generated_per_element) | |
Constructor. More... | |
SamplePointContainer () | |
Broken default constructor; needed for broken copy constructors. Don't call. It will die. More... | |
SamplePointContainer (const SamplePointContainer &data)=delete | |
Broken copy constructor. More... | |
void | operator= (const SamplePointContainer &)=delete |
Broken assignment operator. More... | |
virtual | ~SamplePointContainer () |
Virtual destructor. More... | |
Mesh * | mesh_pt () const |
Pointer to mesh from whose FiniteElements sample points are created. More... | |
const std::pair< double, double > & | min_and_max_coordinates (const unsigned &i) const |
Pair of doubles for min and maximum coordinates in i-th direction: min (first) and max. (second) coordinates. More... | |
const Vector< std::pair< double, double > > & | min_and_max_coordinates () const |
Vector of pair of doubles for min and maximum coordinates. min (first) and max. (second) coordinates. More... | |
bool | ignore_halo_elements_during_locate_zeta_search () const |
Ignore halo elements? More... | |
bool | use_eulerian_coordinates_during_setup () const |
Use Eulerian coordinates (i.e. interpolated_x) rather than zeta itself (i.e. interpolated_zeta) to identify point. More... | |
unsigned & | nsample_points_generated_per_element () |
"Measure of" number of sample points generated in each element More... | |
double & | max_search_radius () |
Set maximum search radius for locate zeta. This is initialised do DBL_MAX so we brutally search through the entire bin structure, no matter how big it is until we've found the required point (or failed to do so. This can be VERY costly with fine meshes. Here the user takes full responsibility and states that we have no chance in hell to find the required point in a bin whose closest vertex is further than the specified max search radius. More... | |
Static Public Attributes | |
static unsigned | Default_n_bin_1d = 5 |
Default number of bins (in each coordinate direction) (Note: don't move this into a common base class because each derived class has its own value; we'll want far fewer in the refineable version!) More... | |
Static Public Attributes inherited from SamplePointContainer | |
static std::ofstream | Visited_sample_points_file |
File to record sequence of visited sample points in. Used for debugging/ illustration of search procedures. More... | |
static bool | Always_fail_elemental_locate_zeta = false |
Boolean flag to make to make locate zeta fail. Used for debugging/ illustration of search procedures. More... | |
static bool | Use_equally_spaced_interior_sample_points = true |
Use equally spaced sample points? (otherwise vertices are sampled repeatedly. More... | |
static bool | Enable_timing_of_setup = false |
Time setup? More... | |
static double | Percentage_offset = 5.0 |
Offset of sample point container boundaries beyond max/min coords. More... | |
Private Member Functions | |
void | fill_bin_array () |
Fill the bin array with sample points from FiniteElements stored in mesh. More... | |
void | create_sample_points_from_element (FiniteElement *const element_pt, const unsigned &n_element) |
Loop over all sample points in the element specified via the pointer and create a SamplePoint for each. Also specify the index of the element in its mesh. More... | |
Private Attributes | |
Vector< RefineableBin * > | Bin_pt |
Vector of pointers to constituent RefineableBins. More... | |
bool | Bin_array_is_recursive |
Variable which stores if the RefineableBinArray is recursive or not. More... | |
unsigned | Depth |
Variable which stores the Depth value of the bin_array. Useful for debugging and for preventing "infinite" recursion in case if there is a problem. More... | |
unsigned | Max_depth |
Max depth of the hierarchical bin_array. More... | |
unsigned | Max_number_of_sample_point_per_bin |
Maximum number of sample points in bin (before it's subdivided recursively) More... | |
RefineableBinArray * | Root_bin_array_pt |
Pointer to root bin array. More... | |
unsigned | First_sample_point_to_actually_lookup_during_locate_zeta |
When searching through sample points recursively from the top level RefineableBinArray (in deterministic order!) only actually do the locate_zeta calls when when the counter exceeds this value. More... | |
unsigned | Last_sample_point_to_actually_lookup_during_locate_zeta |
When searching through sample points recursively from the top level RefineableBinArray (in deterministic order!) only actually do the locate_zeta calls when when the counter is less than this value. More... | |
unsigned | Multiplier_for_max_sample_point_to_actually_lookup_during_locate_zeta |
Every time we've completed a "spiral", visiting a finite number of sample points in a deterministic order, use this multiplier to increase the max. number of sample points to be visited. Using a multiplier rather than a constant increment increases the amount of (more and more unlikely to yield anything!) work done locally before doing another costly mpi round trip when we're already far from the point we're trying to find. More... | |
unsigned | Initial_last_sample_point_to_actually_lookup_during_locate_zeta |
When searching through sample points recursively from the top level RefineableBinArray (in deterministic order!) only actually do the locate_zeta calls when when the counter exceeds this value. This is the initial value when starting the spiral based search. More... | |
Additional Inherited Members | |
Protected Member Functions inherited from SamplePointContainer | |
void | setup_min_and_max_coordinates () |
Helper function to compute the min and max coordinates for the mesh, in each dimension. More... | |
Protected Attributes inherited from BinArray | |
Vector< unsigned > | Dimensions_of_bin_array |
Number of bins in each coordinate direction. More... | |
Protected Attributes inherited from SamplePointContainer | |
Mesh * | Mesh_pt |
Pointer to mesh from whose FiniteElements sample points are created. More... | |
Vector< std::pair< double, double > > | Min_and_max_coordinates |
Vector of pairs of doubles for min and maximum coordinates. Call: Min_and_max_coordinates[j] gives me the pair of min (first) and max. (second) coordinates in the j-th coordinate direction. More... | |
bool | Use_eulerian_coordinates_during_setup |
Use Eulerian coordinates (i.e. interpolated_x) rather than zeta itself (i.e. interpolated_zeta) to identify point. More... | |
bool | Ignore_halo_elements_during_locate_zeta_search |
Ignore halo elements? More... | |
unsigned | Nsample_points_generated_per_element |
"Measure of" number of sample points generated in each element More... | |
unsigned | Total_number_of_sample_points_visited_during_locate_zeta_from_top_level |
Counter to keep track of how many sample points we've visited during top level call to locate_zeta. More... | |
double | Max_search_radius |
Max radius beyond which we stop searching the bin. Initialised to DBL_MAX so keep going until the point is found or until we've searched every single bin. Overwriting this means we won't search in bins whose closest vertex is at a distance greater than Max_search_radius from the point to be located. More... | |
///////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////
RefineableBinArray class.
Definition at line 520 of file sample_point_container.h.
RefineableBinArray::RefineableBinArray | ( | SamplePointContainerParameters * | sample_point_container_parameters_pt | ) |
Constructor.
/////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////////// RefineableBin array class /////////////////////////////////////////////////////////////////////////// ///////////////////////////////////////////////////////////////////////////
Constructor
Definition at line 841 of file sample_point_container.cc.
References Bin_array_is_recursive, Bin_pt, Default_n_bin_1d, Depth, BinArray::Dimensions_of_bin_array, SamplePointContainer::Enable_timing_of_setup, fill_bin_array(), First_sample_point_to_actually_lookup_during_locate_zeta, i, Initial_last_sample_point_to_actually_lookup_during_locate_zeta, Last_sample_point_to_actually_lookup_during_locate_zeta, Max_depth, Max_number_of_sample_point_per_bin, SamplePointContainer::Mesh_pt, SamplePointContainer::Min_and_max_coordinates, Multiplier_for_max_sample_point_to_actually_lookup_during_locate_zeta, oomph::oomph_info, root_bin_array_pt(), Root_bin_array_pt, SamplePointContainer::setup_min_and_max_coordinates(), oomph::TimingHelpers::timer(), total_number_of_sample_points_computed_recursively(), and SamplePointContainer::Total_number_of_sample_points_visited_during_locate_zeta_from_top_level.
|
delete |
Broken copy constructor.
|
inline |
|
inline |
Add specified SamplePoint to RefineableBinArray.
Definition at line 612 of file sample_point_container.h.
References Bin_pt, and BinArray::coords_to_bin_index().
Referenced by fill_bin_array().
|
inline |
Is the BinArray recursive?
Definition at line 659 of file sample_point_container.h.
References Bin_array_is_recursive.
|
inline |
Pointer to i-th bin; can be null if bin is empty.
Definition at line 553 of file sample_point_container.h.
|
private |
Loop over all sample points in the element specified via the pointer and create a SamplePoint for each. Also specify the index of the element in its mesh.
|
inline |
Depth of the hierarchical bin_array.
Definition at line 641 of file sample_point_container.h.
References Depth.
|
private |
Fill the bin array with sample points from FiniteElements stored in mesh.
For all the sample points we have to create ...
Definition at line 1882 of file sample_point_container.cc.
References Bin_pt, BinArray::coords_to_bin_index(), e, i, SamplePointContainer::Mesh_pt, SamplePointContainer::Min_and_max_coordinates, BinArray::ndim_zeta(), SamplePointContainer::Nsample_points_generated_per_element, s, SamplePointContainer::Use_equally_spaced_interior_sample_points, and SamplePointContainer::Use_eulerian_coordinates_during_setup.
Referenced by RefineableBinArray().
|
inline |
Fill the bin array with specified SamplePoints.
Definition at line 574 of file sample_point_container.h.
References add_sample_point(), e, i, SamplePointContainer::Mesh_pt, BinArray::ndim_zeta(), SamplePointContainer::Nsample_points_generated_per_element, s, SamplePointContainer::Use_equally_spaced_interior_sample_points, and SamplePointContainer::Use_eulerian_coordinates_during_setup.
|
inline |
When searching through sample points recursively from the top level RefineableBinArray (in deterministic order!) only actually do the locate_zeta calls when when the counter exceeds this value.
Definition at line 712 of file sample_point_container.h.
References First_sample_point_to_actually_lookup_during_locate_zeta.
Referenced by oomph::Multi_domain_functions::aux_setup_multi_domain_interaction().
void RefineableBinArray::get_bin_boundaries | ( | const unsigned & | bin_index, |
Vector< std::pair< double, double >> & | min_and_max_coordinates | ||
) |
Boundaries of specified bin in each coordinate direction. *.first = min; *.second = max.
Definition at line 995 of file sample_point_container.cc.
References BinArray::dimension_of_bin_array(), BinArray::Dimensions_of_bin_array, SamplePointContainer::Min_and_max_coordinates, and BinArray::ndim_zeta().
Referenced by output_neighbouring_bins().
|
inline |
When searching through sample points recursively from the top level RefineableBinArray (in deterministic order!) only actually do the locate_zeta calls when when the counter exceeds this value. This is the initial value when starting the spiral based search.
Definition at line 741 of file sample_point_container.h.
References Initial_last_sample_point_to_actually_lookup_during_locate_zeta.
Referenced by oomph::Multi_domain_functions::aux_setup_multi_domain_interaction().
|
inline |
When searching through sample points recursively from the top level RefineableBinArray (in deterministic order!) only actually do the locate_zeta calls when when the counter is less than this value.
Definition at line 720 of file sample_point_container.h.
References Last_sample_point_to_actually_lookup_during_locate_zeta.
Referenced by oomph::Multi_domain_functions::aux_setup_multi_domain_interaction().
|
virtual |
Find sub-GeomObject (finite element) and the local coordinate s within it that contains point with global coordinate zeta. sub_geom_object_pt=0 if point can't be found.
Implements SamplePointContainer.
Definition at line 2002 of file sample_point_container.cc.
References SamplePointContainer::Always_fail_elemental_locate_zeta, Bin_pt, BinArray::coords_to_bin_index(), BinArray::coords_to_vectorial_bin_index(), Depth, BinArray::Dimensions_of_bin_array, BinArray::get_neighbouring_bins_helper(), i, SamplePointContainer::max_search_radius(), SamplePointContainer::Min_and_max_coordinates, nbin(), BinArray::ndim_zeta(), s, total_number_of_sample_points_computed_recursively(), and SamplePointContainer::Total_number_of_sample_points_visited_during_locate_zeta_from_top_level.
|
inline |
Max depth of the hierarchical bin_array.
Definition at line 653 of file sample_point_container.h.
References Max_depth.
|
inline |
Max depth of the hierarchical bin_array; const version.
Definition at line 647 of file sample_point_container.h.
References Max_depth.
|
inline |
Maximum number of sample points in bin (before its subdivided recursively)
Definition at line 666 of file sample_point_container.h.
References Max_number_of_sample_point_per_bin.
|
inline |
Every time we've completed a "spiral", visiting a finite number of sample points in a deterministic order, use this multiplier to increase the max. number of sample points to be visited. Using a multiplier rather than a constant increment increases the amount of (more and more unlikely to yield anything!) work done locally before doing another costly mpi round trip when we're already far from the point we're trying to find.
Definition at line 732 of file sample_point_container.h.
References Multiplier_for_max_sample_point_to_actually_lookup_during_locate_zeta.
Referenced by oomph::Multi_domain_functions::aux_setup_multi_domain_interaction().
|
inlinevirtual |
Number of bins (not taking recursion into account)
Implements BinArray.
Definition at line 559 of file sample_point_container.h.
References Bin_pt.
Referenced by locate_zeta(), and total_number_of_sample_points_computed_recursively().
|
delete |
Broken assignment operator.
|
virtual |
Output bin vertices (allowing display of bins as zones).
Output neighbouring bins up to given "radius" of the specified bin.
Implements BinArray.
Definition at line 1027 of file sample_point_container.cc.
|
inlinevirtual |
Output bins.
Loop over bins
Implements BinArray.
Definition at line 672 of file sample_point_container.h.
void RefineableBinArray::output_neighbouring_bins | ( | const unsigned & | bin_index, |
const unsigned & | radius, | ||
std::ofstream & | outfile | ||
) |
Output neighbouring bins at given "radius" of the specified bin.
Output neighbouring bins up to given "radius" of the specified bin.
Definition at line 1044 of file sample_point_container.cc.
References get_bin_boundaries(), BinArray::get_neighbouring_bins_helper(), i, SamplePointContainer::min_and_max_coordinates(), SamplePointContainer::Min_and_max_coordinates, BinArray::ndim_zeta(), and oomph::oomph_info.
|
inline |
Root bin array.
Definition at line 547 of file sample_point_container.h.
References Root_bin_array_pt.
Referenced by RefineableBinArray().
|
virtual |
Compute total number of sample points recursively.
Implements SamplePointContainer.
Definition at line 1819 of file sample_point_container.cc.
References Bin_pt, i, and nbin().
Referenced by oomph::Multi_domain_functions::aux_setup_multi_domain_interaction(), locate_zeta(), and RefineableBinArray().
|
inlinevirtual |
Counter to keep track of how many sample points we've visited during top level call to locate_zeta.
Reimplemented from SamplePointContainer.
Definition at line 696 of file sample_point_container.h.
References Depth, Root_bin_array_pt, SamplePointContainer::Total_number_of_sample_points_visited_during_locate_zeta_from_top_level, and total_number_of_sample_points_visited_during_locate_zeta_from_top_level().
Referenced by total_number_of_sample_points_visited_during_locate_zeta_from_top_level().
|
private |
Variable which stores if the RefineableBinArray is recursive or not.
Definition at line 763 of file sample_point_container.h.
Referenced by bin_array_is_recursive(), and RefineableBinArray().
|
private |
Vector of pointers to constituent RefineableBins.
Definition at line 759 of file sample_point_container.h.
Referenced by add_sample_point(), bin_pt(), fill_bin_array(), locate_zeta(), nbin(), output_bin_vertices(), output_bins(), RefineableBinArray(), total_number_of_sample_points_computed_recursively(), and ~RefineableBinArray().
|
static |
Default number of bins (in each coordinate direction) (Note: don't move this into a common base class because each derived class has its own value; we'll want far fewer in the refineable version!)
Default number of bins (in each coordinate direction)
Definition at line 568 of file sample_point_container.h.
Referenced by RefineableBinArray().
|
private |
Variable which stores the Depth value of the bin_array. Useful for debugging and for preventing "infinite" recursion in case if there is a problem.
Definition at line 768 of file sample_point_container.h.
Referenced by depth(), locate_zeta(), RefineableBinArray(), and total_number_of_sample_points_visited_during_locate_zeta_from_top_level().
|
private |
When searching through sample points recursively from the top level RefineableBinArray (in deterministic order!) only actually do the locate_zeta calls when when the counter exceeds this value.
Definition at line 785 of file sample_point_container.h.
Referenced by first_sample_point_to_actually_lookup_during_locate_zeta(), and RefineableBinArray().
|
private |
When searching through sample points recursively from the top level RefineableBinArray (in deterministic order!) only actually do the locate_zeta calls when when the counter exceeds this value. This is the initial value when starting the spiral based search.
Definition at line 806 of file sample_point_container.h.
Referenced by initial_last_sample_point_to_actually_lookup_during_locate_zeta(), and RefineableBinArray().
|
private |
When searching through sample points recursively from the top level RefineableBinArray (in deterministic order!) only actually do the locate_zeta calls when when the counter is less than this value.
Definition at line 790 of file sample_point_container.h.
Referenced by last_sample_point_to_actually_lookup_during_locate_zeta(), and RefineableBinArray().
|
private |
Max depth of the hierarchical bin_array.
Definition at line 771 of file sample_point_container.h.
Referenced by max_depth(), and RefineableBinArray().
|
private |
Maximum number of sample points in bin (before it's subdivided recursively)
Definition at line 775 of file sample_point_container.h.
Referenced by max_number_of_sample_point_per_bin(), and RefineableBinArray().
|
private |
Every time we've completed a "spiral", visiting a finite number of sample points in a deterministic order, use this multiplier to increase the max. number of sample points to be visited. Using a multiplier rather than a constant increment increases the amount of (more and more unlikely to yield anything!) work done locally before doing another costly mpi round trip when we're already far from the point we're trying to find.
Definition at line 800 of file sample_point_container.h.
Referenced by multiplier_for_max_sample_point_to_actually_lookup_during_locate_zeta(), and RefineableBinArray().
|
private |
Pointer to root bin array.
Definition at line 778 of file sample_point_container.h.
Referenced by RefineableBinArray(), root_bin_array_pt(), and total_number_of_sample_points_visited_during_locate_zeta_from_top_level().