www.mooseframework.org
Classes | Public Types | Public Member Functions | Protected Member Functions | Static Protected Member Functions | Protected Attributes | List of all members
FeatureFloodCount Class Reference

This object will mark nodes or elements of continuous regions all with a unique number for the purpose of counting or "coloring" unique regions in a solution. More...

#include <FeatureFloodCount.h>

Inheritance diagram for FeatureFloodCount:
[legend]

Classes

class  FeatureData
 

Public Types

enum  FieldType {
  FieldType::UNIQUE_REGION, FieldType::VARIABLE_COLORING, FieldType::GHOSTED_ENTITIES, FieldType::HALOS,
  FieldType::CENTROID, FieldType::ACTIVE_BOUNDS
}
 
enum  Status : unsigned char { Status::CLEAR = 0x0, Status::MARKED = 0x1, Status::DIRTY = 0x2, Status::INACTIVE = 0x4 }
 This enumeration is used to indicate status of the grains in the _unique_grains data structure. More...
 

Public Member Functions

 FeatureFloodCount (const InputParameters &parameters)
 
 ~FeatureFloodCount ()
 
virtual void initialSetup () override
 
virtual void meshChanged () override
 
virtual void initialize () override
 
virtual void execute () override
 
virtual void finalize () override
 
virtual Real getValue () override
 
virtual std::size_t getTotalFeatureCount () const
 Returns the total feature count (active and inactive ids, useful for sizing vectors) More...
 
virtual bool doesFeatureIntersectBoundary (unsigned int feature_id) const
 Returns a Boolean indicating whether this feature intersects any boundary. More...
 
virtual const std::vector< unsigned int > & getVarToFeatureVector (dof_id_type elem_id) const
 Returns a list of active unique feature ids for a particular element. More...
 
virtual unsigned int getFeatureVar (unsigned int feature_id) const
 Returns the variable representing the passed in feature. More...
 
std::size_t numCoupledVars () const
 Returns the number of coupled varaibles. More...
 
const std::vector< MooseVariable * > & getCoupledVars () const
 Returns a const vector to the coupled variable pointers. More...
 
virtual Real getEntityValue (dof_id_type entity_id, FieldType field_type, std::size_t var_index=0) const
 
bool isElemental () const
 
const std::vector< FeatureData > & getFeatures () const
 Return a constant reference to the vector of all discovered features. More...
 

Static Public Attributes

static const std::size_t invalid_size_t = std::numeric_limits<std::size_t>::max()
 
static const unsigned int invalid_id = std::numeric_limits<unsigned int>::max()
 

Protected Member Functions

virtual void updateFieldInfo ()
 This method is used to populate any of the data structures used for storing field data (nodal or elemental). More...
 
bool flood (const DofObject *dof_object, std::size_t current_index, FeatureData *feature)
 This method will "mark" all entities on neighboring elements that are above the supplied threshold. More...
 
virtual Real getThreshold (std::size_t current_index) const
 Return the starting comparison threshold to use when inspecting an entity during the flood stage. More...
 
virtual Real getConnectingThreshold (std::size_t current_index) const
 Return the "connecting" comparison threshold to use when inspecting an entity during the flood stage. More...
 
bool compareValueWithThreshold (Real entity_value, Real threshold) const
 This method is used to determine whether the current entity value is part of a feature or not. More...
 
virtual bool isNewFeatureOrConnectedRegion (const DofObject *dof_object, std::size_t &current_index, FeatureData *&feature, Status &status, unsigned int &new_id)
 Method called during the recursive flood routine that should return whether or not the current entity is part of the current feature (if one is being explored), or if it's the start of a new feature. More...
 
void expandPointHalos ()
 This method takes all of the partial features and expands the local, ghosted, and halo sets around those regions to account for the diffuse interface. More...
 
void expandEdgeHalos (unsigned int num_layers_to_expand)
 This method expands the existing halo set by some width determined by the passed in value. More...
 
template<typename T >
void visitNeighborsHelper (const T *curr_entity, std::vector< const T * > neighbor_entities, std::size_t current_index, FeatureData *feature, bool expand_halos_only, bool topological_neighbor, bool disjoint_only)
 The actual logic for visiting neighbors is abstracted out here. More...
 
void prepareDataForTransfer ()
 This routine uses the local flooded data to build up the local feature data structures (_feature_sets). More...
 
void mergeSets ()
 This routine is called on the master rank only and stitches together the partial feature pieces seen on any processor. More...
 
virtual bool areFeaturesMergeable (const FeatureData &f1, const FeatureData &f2) const
 Method for determining whether two features are mergeable. More...
 
void communicateAndMerge ()
 This routine handles all of the serialization, communication and deserialization of the data structures containing FeatureData objects. More...
 
void sortAndLabel ()
 Sort and assign ids to features based on their position in the container after sorting. More...
 
void scatterAndUpdateRanks ()
 Calls buildLocalToGlobalIndices to build the individual local to global indicies for each rank and scatters that information to all ranks. More...
 
virtual void buildLocalToGlobalIndices (std::vector< std::size_t > &local_to_global_all, std::vector< int > &counts) const
 This routine populates a stacked vector of local to global indices per rank and the associated count vector for scattering the vector to the ranks. More...
 
void buildFeatureIdToLocalIndices (unsigned int max_id)
 This method builds a lookup map for retrieving the right local feature (by index) given a global index or id. More...
 
virtual void clearDataStructures ()
 Helper routine for clearing up data structures during initialize and prior to parallel communication. More...
 
void appendPeriodicNeighborNodes (FeatureData &feature) const
 This routine adds the periodic node information to our data structure prior to packing the data this makes those periodic neighbors appear much like ghosted nodes in a multiprocessor setting. More...
 
void updateRegionOffsets ()
 This routine updates the _region_offsets variable which is useful for quickly determining the proper global number for a feature when using multimap mode. More...
 
void visitNodalNeighbors (const Node *node, std::size_t current_index, FeatureData *feature, bool expand_halos_only)
 These two routines are utility routines used by the flood routine and by derived classes for visiting neighbors. More...
 
void visitElementalNeighbors (const Elem *elem, std::size_t current_index, FeatureData *feature, bool expand_halos_only, bool disjoint_only)
 
void serialize (std::string &serialized_buffer)
 These routines packs/unpack the _feature_map data into a structure suitable for parallel communication operations. More...
 
void deserialize (std::vector< std::string > &serialized_buffers)
 This routine takes the vector of byte buffers (one for each processor), deserializes them into a series of FeatureSet objects, and appends them to the _feature_sets data structure. More...
 

Static Protected Member Functions

template<class InputIterator >
static bool setsIntersect (InputIterator first1, InputIterator last1, InputIterator first2, InputIterator last2)
 This method detects whether two sets intersect without building a result set. More...
 

Protected Attributes

std::vector< MooseVariable * > _vars
 The vector of coupled in variables. More...
 
const Real _threshold
 The threshold above (or below) where an entity may begin a new region (feature) More...
 
Real _step_threshold
 
const Real _connecting_threshold
 The threshold above (or below) which neighboring entities are flooded (where regions can be extended but not started) More...
 
Real _step_connecting_threshold
 
MooseMesh & _mesh
 A reference to the mesh. More...
 
unsigned long _var_number
 This variable is used to build the periodic node map. More...
 
const bool _single_map_mode
 This variable is used to indicate whether or not multiple maps are used during flooding. More...
 
const bool _condense_map_info
 
const bool _global_numbering
 This variable is used to indicate whether or not we identify features with unique numbers on multiple maps. More...
 
const bool _var_index_mode
 This variable is used to indicate whether the maps will contain unique region information or just the variable numbers owning those regions. More...
 
const bool _compute_halo_maps
 Indicates whether or not to communicate halo map information with all ranks. More...
 
const bool _compute_var_to_feature_map
 Indicates whether or not the var to feature map is populated. More...
 
const bool _use_less_than_threshold_comparison
 Use less-than when comparing values against the threshold value. More...
 
const std::size_t _n_vars
 
const std::size_t _maps_size
 Convenience variable holding the size of all the datastructures size by the number of maps. More...
 
const processor_id_type _n_procs
 Convenience variable holding the number of processors in this simulation. More...
 
std::vector< std::set< dof_id_type > > _entities_visited
 This variable keeps track of which nodes have been visited during execution. More...
 
std::vector< std::map< dof_id_type, int > > _var_index_maps
 This map keeps track of which variables own which nodes. More...
 
std::vector< std::vector< const Elem * > > _nodes_to_elem_map
 The data structure used to find neighboring elements give a node ID. More...
 
std::vector< unsigned int > _feature_counts_per_map
 The number of features seen by this object per map. More...
 
unsigned int _feature_count
 The number of features seen by this object (same as summing _feature_counts_per_map) More...
 
std::vector< std::list< FeatureData > > _partial_feature_sets
 The data structure used to hold partial and communicated feature data. More...
 
std::vector< FeatureData_feature_sets
 The data structure used to hold the globally unique features. More...
 
std::vector< std::map< dof_id_type, int > > _feature_maps
 The feature maps contain the raw flooded node information and eventually the unique grain numbers. More...
 
std::vector< std::size_t > _local_to_global_feature_map
 The vector recording the local to global feature indices. More...
 
std::vector< std::size_t > _feature_id_to_local_index
 The vector recording the grain_id to local index (several indices will contain invalid_size_t) More...
 
PeriodicBoundaries * _pbs
 A pointer to the periodic boundary constraints object. More...
 
std::unique_ptr< PointLocatorBase > _point_locator
 
const PostprocessorValue & _element_average_value
 Average value of the domain which can optionally be used to find features in a field. More...
 
std::map< dof_id_type, int > _ghosted_entity_ids
 The map for holding reconstructed ghosted element information. More...
 
std::vector< std::map< dof_id_type, int > > _halo_ids
 The data structure for looking up halos around features. More...
 
std::multimap< dof_id_type, dof_id_type > _periodic_node_map
 The data structure which is a list of nodes that are constrained to other nodes based on the imposed periodic boundary conditions. More...
 
std::set< dof_id_type > _all_boundary_entity_ids
 The set of entities on the boundary of the domain used for determining if features intersect any boundary. More...
 
std::map< dof_id_type, std::vector< unsigned int > > _entity_var_to_features
 
std::vector< unsigned int > _empty_var_to_features
 
bool _is_elemental
 Determines if the flood counter is elements or not (nodes) More...
 
bool _is_master
 Convenience variable for testing master rank. More...
 

Detailed Description

This object will mark nodes or elements of continuous regions all with a unique number for the purpose of counting or "coloring" unique regions in a solution.

It is designed to work with either a single variable, or multiple variables.

Note: When inspecting multiple variables, those variables must not have regions of interest that overlap or they will not be correctly colored.

Definition at line 43 of file FeatureFloodCount.h.

Member Enumeration Documentation

Enumerator
UNIQUE_REGION 
VARIABLE_COLORING 
GHOSTED_ENTITIES 
HALOS 
CENTROID 
ACTIVE_BOUNDS 

Definition at line 88 of file FeatureFloodCount.h.

89  {
90  UNIQUE_REGION,
91  VARIABLE_COLORING,
92  GHOSTED_ENTITIES,
93  HALOS,
94  CENTROID,
95  ACTIVE_BOUNDS,
96  };
enum FeatureFloodCount::Status : unsigned char
strong

This enumeration is used to indicate status of the grains in the _unique_grains data structure.

Enumerator
CLEAR 
MARKED 
DIRTY 
INACTIVE 

Definition at line 105 of file FeatureFloodCount.h.

105  : unsigned char
106  {
107  CLEAR = 0x0,
108  MARKED = 0x1,
109  DIRTY = 0x2,
110  INACTIVE = 0x4
111  };

Constructor & Destructor Documentation

FeatureFloodCount::FeatureFloodCount ( const InputParameters &  parameters)

Definition at line 156 of file FeatureFloodCount.C.

157  : GeneralPostprocessor(parameters),
158  Coupleable(this, false),
159  MooseVariableDependencyInterface(),
160  ZeroInterface(parameters),
161  _vars(getCoupledMooseVars()),
162  _threshold(getParam<Real>("threshold")),
163  _connecting_threshold(isParamValid("connecting_threshold")
164  ? getParam<Real>("connecting_threshold")
165  : getParam<Real>("threshold")),
166  _mesh(_subproblem.mesh()),
167  _var_number(_vars[0]->number()),
168  _single_map_mode(getParam<bool>("use_single_map")),
169  _condense_map_info(getParam<bool>("condense_map_info")),
170  _global_numbering(getParam<bool>("use_global_numbering")),
171  _var_index_mode(getParam<bool>("enable_var_coloring")),
172  _compute_halo_maps(getParam<bool>("compute_halo_maps")),
173  _compute_var_to_feature_map(getParam<bool>("compute_var_to_feature_map")),
174  _use_less_than_threshold_comparison(getParam<bool>("use_less_than_threshold_comparison")),
175  _n_vars(_vars.size()),
176  _maps_size(_single_map_mode ? 1 : _vars.size()),
177  _n_procs(_app.n_processors()),
179  _feature_count(0),
182  _pbs(nullptr),
183  _element_average_value(parameters.isParamValid("elem_avg_value")
184  ? getPostprocessorValue("elem_avg_value")
185  : _real_zero),
187  _is_elemental(getParam<MooseEnum>("flood_entity_type") == "ELEMENTAL"),
188  _is_master(processor_id() == 0)
189 {
190  if (_var_index_mode)
191  _var_index_maps.resize(_maps_size);
192 
193  addMooseVariableDependency(_vars);
194 }
const PostprocessorValue & _element_average_value
Average value of the domain which can optionally be used to find features in a field.
const std::size_t _n_vars
const bool _condense_map_info
std::vector< std::map< dof_id_type, int > > _halo_ids
The data structure for looking up halos around features.
std::vector< MooseVariable * > _vars
The vector of coupled in variables.
const Real _connecting_threshold
The threshold above (or below) which neighboring entities are flooded (where regions can be extended ...
bool _is_elemental
Determines if the flood counter is elements or not (nodes)
bool _is_master
Convenience variable for testing master rank.
unsigned long _var_number
This variable is used to build the periodic node map.
const std::size_t _maps_size
Convenience variable holding the size of all the datastructures size by the number of maps...
const bool _use_less_than_threshold_comparison
Use less-than when comparing values against the threshold value.
const bool _global_numbering
This variable is used to indicate whether or not we identify features with unique numbers on multiple...
std::vector< std::map< dof_id_type, int > > _feature_maps
The feature maps contain the raw flooded node information and eventually the unique grain numbers...
const bool _single_map_mode
This variable is used to indicate whether or not multiple maps are used during flooding.
unsigned int _feature_count
The number of features seen by this object (same as summing _feature_counts_per_map) ...
const Real _threshold
The threshold above (or below) where an entity may begin a new region (feature)
std::vector< std::list< FeatureData > > _partial_feature_sets
The data structure used to hold partial and communicated feature data.
const bool _compute_halo_maps
Indicates whether or not to communicate halo map information with all ranks.
PeriodicBoundaries * _pbs
A pointer to the periodic boundary constraints object.
const processor_id_type _n_procs
Convenience variable holding the number of processors in this simulation.
std::vector< unsigned int > _feature_counts_per_map
The number of features seen by this object per map.
MooseMesh & _mesh
A reference to the mesh.
const bool _compute_var_to_feature_map
Indicates whether or not the var to feature map is populated.
const bool _var_index_mode
This variable is used to indicate whether the maps will contain unique region information or just the...
std::vector< std::map< dof_id_type, int > > _var_index_maps
This map keeps track of which variables own which nodes.
FeatureFloodCount::~FeatureFloodCount ( )

Definition at line 196 of file FeatureFloodCount.C.

196 {}

Member Function Documentation

void FeatureFloodCount::appendPeriodicNeighborNodes ( FeatureData feature) const
protected

This routine adds the periodic node information to our data structure prior to packing the data this makes those periodic neighbors appear much like ghosted nodes in a multiprocessor setting.

Definition at line 1382 of file FeatureFloodCount.C.

Referenced by getFeatures(), and prepareDataForTransfer().

1383 {
1384  if (_is_elemental)
1385  {
1386  for (auto entity : feature._local_ids)
1387  {
1388  Elem * elem = _mesh.elemPtr(entity);
1389 
1390  for (auto node_n = decltype(elem->n_nodes())(0); node_n < elem->n_nodes(); ++node_n)
1391  {
1392  auto iters = _periodic_node_map.equal_range(elem->node(node_n));
1393 
1394  for (auto it = iters.first; it != iters.second; ++it)
1395  {
1396  feature._periodic_nodes.insert(it->first);
1397  feature._periodic_nodes.insert(it->second);
1398  }
1399  }
1400  }
1401  }
1402  else
1403  {
1404  for (auto entity : feature._local_ids)
1405  {
1406  auto iters = _periodic_node_map.equal_range(entity);
1407 
1408  for (auto it = iters.first; it != iters.second; ++it)
1409  {
1410  feature._periodic_nodes.insert(it->first);
1411  feature._periodic_nodes.insert(it->second);
1412  }
1413  }
1414  }
1415 }
std::multimap< dof_id_type, dof_id_type > _periodic_node_map
The data structure which is a list of nodes that are constrained to other nodes based on the imposed ...
bool _is_elemental
Determines if the flood counter is elements or not (nodes)
MooseMesh & _mesh
A reference to the mesh.
bool FeatureFloodCount::areFeaturesMergeable ( const FeatureData f1,
const FeatureData f2 
) const
protectedvirtual

Method for determining whether two features are mergeable.

This routine exists because derived classes may need to override this function rather than use the mergeable method in the FeatureData object.

Reimplemented in PolycrystalUserObjectBase.

Definition at line 933 of file FeatureFloodCount.C.

Referenced by getFeatures(), and mergeSets().

934 {
935  return f1.mergeable(f2);
936 }
void FeatureFloodCount::buildFeatureIdToLocalIndices ( unsigned int  max_id)
protected

This method builds a lookup map for retrieving the right local feature (by index) given a global index or id.

max_id is passed to size the vector properly and may or may not be a globally consistent number. The assumption is that any id that is later queried from this object that is higher simply doesn't exist on the local processor.

Definition at line 453 of file FeatureFloodCount.C.

Referenced by GrainTracker::assignGrains(), getFeatures(), scatterAndUpdateRanks(), and GrainTracker::trackGrains().

454 {
455  _feature_id_to_local_index.assign(max_id + 1, invalid_size_t);
456  for (auto feature_index = beginIndex(_feature_sets); feature_index < _feature_sets.size();
457  ++feature_index)
458  {
459  if (_feature_sets[feature_index]._status != Status::INACTIVE)
460  {
461  mooseAssert(_feature_sets[feature_index]._id <= max_id,
462  "Feature ID out of range(" << _feature_sets[feature_index]._id << ')');
463  _feature_id_to_local_index[_feature_sets[feature_index]._id] = feature_index;
464  }
465  }
466 }
static const std::size_t invalid_size_t
std::vector< FeatureData > _feature_sets
The data structure used to hold the globally unique features.
std::vector< std::size_t > _feature_id_to_local_index
The vector recording the grain_id to local index (several indices will contain invalid_size_t) ...
void FeatureFloodCount::buildLocalToGlobalIndices ( std::vector< std::size_t > &  local_to_global_all,
std::vector< int > &  counts 
) const
protectedvirtual

This routine populates a stacked vector of local to global indices per rank and the associated count vector for scattering the vector to the ranks.

The individual vectors can be different sizes. The ith vector will be distributed to the ith processor including the master rank. e.g. [ ... n_0 ] [ ... n_1 ] ... [ ... n_m ]

It is intended to be overridden in derived classes.

Definition at line 407 of file FeatureFloodCount.C.

Referenced by getFeatures(), and scatterAndUpdateRanks().

409 {
410  mooseAssert(_is_master, "This method must only be called on the root processor");
411 
412  counts.assign(_n_procs, 0);
413  // Now size the individual counts vectors based on the largest index seen per processor
414  for (const auto & feature : _feature_sets)
415  for (const auto & local_index_pair : feature._orig_ids)
416  {
417  // local_index_pair.first = ranks, local_index_pair.second = local_index
418  mooseAssert(local_index_pair.first < _n_procs, "Processor ID is out of range");
419  if (local_index_pair.second >= static_cast<std::size_t>(counts[local_index_pair.first]))
420  counts[local_index_pair.first] = local_index_pair.second + 1;
421  }
422 
423  // Build the offsets vector
424  unsigned int globalsize = 0;
425  std::vector<int> offsets(_n_procs); // Type is signed for use with the MPI API
426  for (auto i = beginIndex(offsets); i < offsets.size(); ++i)
427  {
428  offsets[i] = globalsize;
429  globalsize += counts[i];
430  }
431 
432  // Finally populate the master vector
433  local_to_global_all.resize(globalsize, FeatureFloodCount::invalid_size_t);
434  for (const auto & feature : _feature_sets)
435  {
436  // Get the local indices from the feature and build a map
437  for (const auto & local_index_pair : feature._orig_ids)
438  {
439  auto rank = local_index_pair.first;
440  mooseAssert(rank < _n_procs, rank << ", " << _n_procs);
441 
442  auto local_index = local_index_pair.second;
443  auto stacked_local_index = offsets[rank] + local_index;
444 
445  mooseAssert(stacked_local_index < globalsize,
446  "Global index: " << stacked_local_index << " is out of range");
447  local_to_global_all[stacked_local_index] = feature._id;
448  }
449  }
450 }
static const std::size_t invalid_size_t
bool _is_master
Convenience variable for testing master rank.
std::vector< FeatureData > _feature_sets
The data structure used to hold the globally unique features.
const processor_id_type _n_procs
Convenience variable holding the number of processors in this simulation.
void FeatureFloodCount::clearDataStructures ( )
protectedvirtual

Helper routine for clearing up data structures during initialize and prior to parallel communication.

Definition at line 252 of file FeatureFloodCount.C.

Referenced by communicateAndMerge(), and getFeatures().

253 {
254 }
void FeatureFloodCount::communicateAndMerge ( )
protected

This routine handles all of the serialization, communication and deserialization of the data structures containing FeatureData objects.

The libMesh packed range routines handle the communication of the individual string buffers. Here we need to create a container to hold our type to serialize. It'll always be size one because we are sending a single byte stream of all the data to other processors. The stream need not be the same size on all processors.

Additionally we need to create a different container to hold the received byte buffers. The container type need not match the send container type. However, We do know the number of incoming buffers (num processors) so we'll go ahead and use a vector.

Send the data from all processors to the root to create a complete global feature map.

Definition at line 308 of file FeatureFloodCount.C.

Referenced by GrainTracker::finalize(), finalize(), and getFeatures().

309 {
310  // First we need to transform the raw data into a usable data structure
312 
320  std::vector<std::string> send_buffers(1);
321 
328  std::vector<std::string> recv_buffers;
329  if (_is_master)
330  recv_buffers.reserve(_app.n_processors());
331 
332  serialize(send_buffers[0]);
333 
334  // Free up as much memory as possible here before we do global communication
336 
341  _communicator.gather_packed_range(0,
342  (void *)(nullptr),
343  send_buffers.begin(),
344  send_buffers.end(),
345  std::back_inserter(recv_buffers));
346 
347  if (_is_master)
348  {
349  // The root process now needs to deserialize and merge all of the data
350  deserialize(recv_buffers);
351  recv_buffers.clear();
352 
353  mergeSets();
354  }
355 
356  // Make sure that feature count is communicated to all ranks
357  _communicator.broadcast(_feature_count);
358 }
void serialize(std::string &serialized_buffer)
These routines packs/unpack the _feature_map data into a structure suitable for parallel communicatio...
void deserialize(std::vector< std::string > &serialized_buffers)
This routine takes the vector of byte buffers (one for each processor), deserializes them into a seri...
bool _is_master
Convenience variable for testing master rank.
virtual void clearDataStructures()
Helper routine for clearing up data structures during initialize and prior to parallel communication...
unsigned int _feature_count
The number of features seen by this object (same as summing _feature_counts_per_map) ...
void mergeSets()
This routine is called on the master rank only and stitches together the partial feature pieces seen ...
void prepareDataForTransfer()
This routine uses the local flooded data to build up the local feature data structures (_feature_sets...
bool FeatureFloodCount::compareValueWithThreshold ( Real  entity_value,
Real  threshold 
) const
protected

This method is used to determine whether the current entity value is part of a feature or not.

Comparisons can either be greater than or less than the threshold which is controlled via input parameter.

Definition at line 1083 of file FeatureFloodCount.C.

Referenced by getFeatures(), and isNewFeatureOrConnectedRegion().

1084 {
1085  return ((_use_less_than_threshold_comparison && (entity_value >= threshold)) ||
1086  (!_use_less_than_threshold_comparison && (entity_value <= threshold)));
1087 }
const bool _use_less_than_threshold_comparison
Use less-than when comparing values against the threshold value.
void FeatureFloodCount::deserialize ( std::vector< std::string > &  serialized_buffers)
protected

This routine takes the vector of byte buffers (one for each processor), deserializes them into a series of FeatureSet objects, and appends them to the _feature_sets data structure.

Note: It is assumed that local processor information may already be stored in the _feature_sets data structure so it is not cleared before insertion.

We should already have the local processor data in the features data structure. Don't unpack the local buffer again.

Definition at line 801 of file FeatureFloodCount.C.

Referenced by communicateAndMerge(), and getFeatures().

802 {
803  // The input string stream used for deserialization
804  std::istringstream iss;
805 
806  mooseAssert(serialized_buffers.size() == _app.n_processors(),
807  "Unexpected size of serialized_buffers: " << serialized_buffers.size());
808  auto rank = processor_id();
809  for (auto proc_id = beginIndex(serialized_buffers); proc_id < serialized_buffers.size();
810  ++proc_id)
811  {
816  if (proc_id == rank)
817  continue;
818 
819  iss.str(serialized_buffers[proc_id]); // populate the stream with a new buffer
820  iss.clear(); // reset the string stream state
821 
822  // Load the communicated data into all of the other processors' slots
823  dataLoad(iss, _partial_feature_sets, this);
824  }
825 }
void dataLoad(std::istream &stream, FeatureFloodCount::FeatureData &feature, void *context)
std::vector< std::list< FeatureData > > _partial_feature_sets
The data structure used to hold partial and communicated feature data.
bool FeatureFloodCount::doesFeatureIntersectBoundary ( unsigned int  feature_id) const
virtual

Returns a Boolean indicating whether this feature intersects any boundary.

Reimplemented in GrainTracker, and FauxGrainTracker.

Definition at line 606 of file FeatureFloodCount.C.

Referenced by FeatureVolumeVectorPostprocessor::execute().

607 {
608  // TODO: This information is not parallel consistent when using FeatureFloodCounter
609 
610  // Some processors don't contain the largest feature id, in that case we just return invalid_id
611  if (feature_id >= _feature_id_to_local_index.size())
612  return false;
613 
614  auto local_index = _feature_id_to_local_index[feature_id];
615 
616  if (local_index != invalid_size_t)
617  {
618  mooseAssert(local_index < _feature_sets.size(), "local_index out of bounds");
619  return _feature_sets[local_index]._intersects_boundary;
620  }
621 
622  return false;
623 }
static const std::size_t invalid_size_t
std::vector< FeatureData > _feature_sets
The data structure used to hold the globally unique features.
std::vector< std::size_t > _feature_id_to_local_index
The vector recording the grain_id to local index (several indices will contain invalid_size_t) ...
void FeatureFloodCount::execute ( )
overridevirtual

Reimplemented in PolycrystalUserObjectBase, FauxGrainTracker, and GrainTracker.

Definition at line 280 of file FeatureFloodCount.C.

Referenced by GrainTracker::execute().

281 {
282  const auto end = _mesh.getMesh().active_local_elements_end();
283  for (auto el = _mesh.getMesh().active_local_elements_begin(); el != end; ++el)
284  {
285  const Elem * current_elem = *el;
286 
287  // Loop over elements or nodes
288  if (_is_elemental)
289  {
290  for (auto var_num = beginIndex(_vars); var_num < _vars.size(); ++var_num)
291  flood(current_elem, var_num, nullptr /* Designates inactive feature */);
292  }
293  else
294  {
295  auto n_nodes = current_elem->n_vertices();
296  for (auto i = decltype(n_nodes)(0); i < n_nodes; ++i)
297  {
298  const Node * current_node = current_elem->get_node(i);
299 
300  for (auto var_num = beginIndex(_vars); var_num < _vars.size(); ++var_num)
301  flood(current_node, var_num, nullptr /* Designates inactive feature */);
302  }
303  }
304  }
305 }
std::vector< MooseVariable * > _vars
The vector of coupled in variables.
bool _is_elemental
Determines if the flood counter is elements or not (nodes)
bool flood(const DofObject *dof_object, std::size_t current_index, FeatureData *feature)
This method will "mark" all entities on neighboring elements that are above the supplied threshold...
MooseMesh & _mesh
A reference to the mesh.
void FeatureFloodCount::expandEdgeHalos ( unsigned int  num_layers_to_expand)
protected

This method expands the existing halo set by some width determined by the passed in value.

This method does NOT mask off any local IDs.

Create a copy of the halo set so that as we insert new ids into the set we don't continue to iterate on those new ids.

We have to handle disjoint halo IDs slightly differently. Once you are disjoint, you can't go back so make sure that we keep placing these IDs in the disjoint set.

Definition at line 1189 of file FeatureFloodCount.C.

Referenced by GrainTracker::finalize(), PolycrystalUserObjectBase::finalize(), and getFeatures().

1190 {
1191  if (num_layers_to_expand == 0)
1192  return;
1193 
1194  for (auto & list_ref : _partial_feature_sets)
1195  {
1196  for (auto & feature : list_ref)
1197  {
1198  for (auto halo_level = decltype(num_layers_to_expand)(0); halo_level < num_layers_to_expand;
1199  ++halo_level)
1200  {
1205  std::set<dof_id_type> orig_halo_ids(feature._halo_ids);
1206  for (auto entity : orig_halo_ids)
1207  {
1208  if (_is_elemental)
1209  visitElementalNeighbors(_mesh.elemPtr(entity),
1210  feature._var_index,
1211  &feature,
1212  /*expand_halos_only =*/true,
1213  /*disjoint_only =*/false);
1214  else
1215  visitNodalNeighbors(_mesh.nodePtr(entity),
1216  feature._var_index,
1217  &feature,
1218  /*expand_halos_only =*/true);
1219  }
1220 
1225  std::set<dof_id_type> disjoint_orig_halo_ids(feature._disjoint_halo_ids);
1226  for (auto entity : disjoint_orig_halo_ids)
1227  {
1228  if (_is_elemental)
1229  visitElementalNeighbors(_mesh.elemPtr(entity),
1230  feature._var_index,
1231  &feature,
1232  /*expand_halos_only =*/true,
1233  /*disjoint_only =*/true);
1234  else
1235  visitNodalNeighbors(_mesh.nodePtr(entity),
1236  feature._var_index,
1237  &feature,
1238  /*expand_halos_only =*/true);
1239  }
1240  }
1241  }
1242  }
1243 }
void visitElementalNeighbors(const Elem *elem, std::size_t current_index, FeatureData *feature, bool expand_halos_only, bool disjoint_only)
void visitNodalNeighbors(const Node *node, std::size_t current_index, FeatureData *feature, bool expand_halos_only)
These two routines are utility routines used by the flood routine and by derived classes for visiting...
bool _is_elemental
Determines if the flood counter is elements or not (nodes)
std::vector< std::list< FeatureData > > _partial_feature_sets
The data structure used to hold partial and communicated feature data.
MooseMesh & _mesh
A reference to the mesh.
void FeatureFloodCount::expandPointHalos ( )
protected

This method takes all of the partial features and expands the local, ghosted, and halo sets around those regions to account for the diffuse interface.

Rather than using any kind of recursion here, we simply expand the region by all "point" neighbors from the actual grain cells since all point neighbors will contain contributions to the region.

To expand the feature element region to the actual flooded region (nodal basis) we need to add in all point neighbors of the current local region for each feature. This is because the elemental variable influence spreads from the elemental data out exactly one element from every mesh point.

Definition at line 1130 of file FeatureFloodCount.C.

Referenced by getFeatures().

1131 {
1132  const auto & node_to_elem_map = _mesh.nodeToActiveSemilocalElemMap();
1133  decltype(FeatureData::_local_ids) expanded_local_ids;
1134  auto my_processor_id = processor_id();
1135 
1142  for (auto & list_ref : _partial_feature_sets)
1143  {
1144  for (auto & feature : list_ref)
1145  {
1146  expanded_local_ids.clear();
1147 
1148  for (auto entity : feature._local_ids)
1149  {
1150  const Elem * elem = _mesh.elemPtr(entity);
1151  mooseAssert(elem, "elem pointer is NULL");
1152 
1153  // Get the nodes on a current element so that we can add in point neighbors
1154  auto n_nodes = elem->n_vertices();
1155  for (auto i = decltype(n_nodes)(0); i < n_nodes; ++i)
1156  {
1157  const Node * current_node = elem->get_node(i);
1158 
1159  auto elem_vector_it = node_to_elem_map.find(current_node->id());
1160  if (elem_vector_it == node_to_elem_map.end())
1161  mooseError("Error in node to elem map");
1162 
1163  const auto & elem_vector = elem_vector_it->second;
1164 
1165  expanded_local_ids.insert(elem_vector.begin(), elem_vector.end());
1166 
1167  // Now see which elements need to go into the ghosted set
1168  for (auto entity : elem_vector)
1169  {
1170  const Elem * neighbor = _mesh.elemPtr(entity);
1171  mooseAssert(neighbor, "neighbor pointer is NULL");
1172 
1173  if (neighbor->processor_id() != my_processor_id)
1174  feature._ghosted_ids.insert(elem->id());
1175  }
1176  }
1177  }
1178 
1179  // Replace the existing local ids with the expanded local ids
1180  feature._local_ids.swap(expanded_local_ids);
1181 
1182  // Copy the expanded local_ids into the halo_ids container
1183  feature._halo_ids = feature._local_ids;
1184  }
1185  }
1186 }
std::set< dof_id_type > _local_ids
Holds the local ids in the interior of a feature.
std::vector< std::list< FeatureData > > _partial_feature_sets
The data structure used to hold partial and communicated feature data.
MooseMesh & _mesh
A reference to the mesh.
void FeatureFloodCount::finalize ( )
overridevirtual

Reimplemented in PolycrystalUserObjectBase, FauxGrainTracker, and GrainTracker.

Definition at line 469 of file FeatureFloodCount.C.

Referenced by PolycrystalUserObjectBase::finalize().

470 {
471  // Gather all information on processor zero and merge
473 
474  // Sort and label the features
475  if (_is_master)
476  sortAndLabel();
477 
478  // Send out the local to global mappings
480 
481  // Populate _feature_maps and _var_index_maps
482  updateFieldInfo();
483 }
virtual void updateFieldInfo()
This method is used to populate any of the data structures used for storing field data (nodal or elem...
void communicateAndMerge()
This routine handles all of the serialization, communication and deserialization of the data structur...
void sortAndLabel()
Sort and assign ids to features based on their position in the container after sorting.
bool _is_master
Convenience variable for testing master rank.
void scatterAndUpdateRanks()
Calls buildLocalToGlobalIndices to build the individual local to global indicies for each rank and sc...
bool FeatureFloodCount::flood ( const DofObject *  dof_object,
std::size_t  current_index,
FeatureData feature 
)
protected

This method will "mark" all entities on neighboring elements that are above the supplied threshold.

If feature is NULL, we are exploring for a new region to mark, otherwise we are in the recursive calls currently marking a region.

Returns
Boolean indicating whether a new feature was found while exploring the current entity.

If we reach this point (i.e. we haven't returned early from this routine), we've found a new mesh entity that's part of a feature. We need to mark the entity as visited at this point (and not before!) to avoid infinite recursion. If you mark the node too early you risk not coloring in a whole feature any time a "connecting threshold" is used since we may have already visited this entity earlier but it was in-between two thresholds.

See if this particular entity cell contributes to the centroid calculation. We only deal with elemental floods and only count it if it's owned by the current processor to avoid skewing the result.

Definition at line 984 of file FeatureFloodCount.C.

Referenced by execute(), PolycrystalUserObjectBase::execute(), getFeatures(), and visitNeighborsHelper().

987 {
988  if (dof_object == nullptr)
989  return false;
990 
991  // Retrieve the id of the current entity
992  auto entity_id = dof_object->id();
993 
994  // Has this entity already been marked? - if so move along
995  if (current_index != invalid_size_t &&
996  _entities_visited[current_index].find(entity_id) != _entities_visited[current_index].end())
997  return false;
998 
999  // See if the current entity either starts a new feature or continues an existing feature
1000  auto new_id = invalid_id; // Writable reference to hold an optional id;
1001  Status status =
1002  Status::INACTIVE; // Status is inactive until we find an entity above the starting threshold
1003  if (!isNewFeatureOrConnectedRegion(dof_object, current_index, feature, status, new_id))
1004  return false;
1005 
1006  mooseAssert(current_index != invalid_size_t, "current_index is invalid");
1007 
1016  _entities_visited[current_index].emplace(entity_id);
1017 
1018  auto map_num = _single_map_mode ? decltype(current_index)(0) : current_index;
1019 
1020  // New Feature (we need to create it and add it to our data structure)
1021  if (!feature)
1022  {
1023  _partial_feature_sets[map_num].emplace_back(
1024  current_index, _feature_count++, processor_id(), status);
1025 
1026  // Get a handle to the feature we will update (always the last feature in the data structure)
1027  feature = &_partial_feature_sets[map_num].back();
1028 
1029  // If new_id is valid, we'll set it in the feature here.
1030  if (new_id != invalid_id)
1031  feature->_id = new_id;
1032  }
1033 
1034  // Insert the current entity into the local ids map
1035  feature->_local_ids.insert(entity_id);
1036 
1042  if (_is_elemental && processor_id() == dof_object->processor_id())
1043  {
1044  const Elem * elem = static_cast<const Elem *>(dof_object);
1045 
1046  // Keep track of how many elements participate in the centroid averaging
1047  feature->_vol_count++;
1048 
1049  // Sum the centroid values for now, we'll average them later
1050  feature->_centroid += elem->centroid();
1051 
1052  // Does the volume intersect the boundary?
1053  if (_all_boundary_entity_ids.find(elem->id()) != _all_boundary_entity_ids.end())
1054  feature->_intersects_boundary = true;
1055  }
1056 
1057  if (_is_elemental)
1058  visitElementalNeighbors(static_cast<const Elem *>(dof_object),
1059  current_index,
1060  feature,
1061  /*expand_halos_only =*/false,
1062  /*disjoint_only =*/false);
1063  else
1064  visitNodalNeighbors(static_cast<const Node *>(dof_object),
1065  current_index,
1066  feature,
1067  /*expand_halos_only =*/false);
1068 
1069  return true;
1070 }
void visitElementalNeighbors(const Elem *elem, std::size_t current_index, FeatureData *feature, bool expand_halos_only, bool disjoint_only)
static const std::size_t invalid_size_t
Status
This enumeration is used to indicate status of the grains in the _unique_grains data structure...
std::vector< std::set< dof_id_type > > _entities_visited
This variable keeps track of which nodes have been visited during execution.
void visitNodalNeighbors(const Node *node, std::size_t current_index, FeatureData *feature, bool expand_halos_only)
These two routines are utility routines used by the flood routine and by derived classes for visiting...
bool _is_elemental
Determines if the flood counter is elements or not (nodes)
std::set< dof_id_type > _all_boundary_entity_ids
The set of entities on the boundary of the domain used for determining if features intersect any boun...
static const unsigned int invalid_id
const bool _single_map_mode
This variable is used to indicate whether or not multiple maps are used during flooding.
unsigned int _feature_count
The number of features seen by this object (same as summing _feature_counts_per_map) ...
virtual bool isNewFeatureOrConnectedRegion(const DofObject *dof_object, std::size_t &current_index, FeatureData *&feature, Status &status, unsigned int &new_id)
Method called during the recursive flood routine that should return whether or not the current entity...
std::vector< std::list< FeatureData > > _partial_feature_sets
The data structure used to hold partial and communicated feature data.
Real FeatureFloodCount::getConnectingThreshold ( std::size_t  current_index) const
protectedvirtual

Return the "connecting" comparison threshold to use when inspecting an entity during the flood stage.

Definition at line 1077 of file FeatureFloodCount.C.

Referenced by getFeatures(), and isNewFeatureOrConnectedRegion().

1078 {
1080 }
const std::vector<MooseVariable *>& FeatureFloodCount::getCoupledVars ( ) const
inline

Returns a const vector to the coupled variable pointers.

Definition at line 86 of file FeatureFloodCount.h.

Referenced by AverageGrainVolume::AverageGrainVolume().

86 { return _vars; }
std::vector< MooseVariable * > _vars
The vector of coupled in variables.
Real FeatureFloodCount::getEntityValue ( dof_id_type  entity_id,
FieldType  field_type,
std::size_t  var_index = 0 
) const
virtual

Reimplemented in GrainTracker, and FauxGrainTracker.

Definition at line 626 of file FeatureFloodCount.C.

Referenced by GrainTracker::getEntityValue(), and FeatureFloodCountAux::precalculateValue().

629 {
630  auto use_default = false;
631  if (var_index == invalid_size_t)
632  {
633  use_default = true;
634  var_index = 0;
635  }
636 
637  mooseAssert(var_index < _maps_size, "Index out of range");
638 
639  switch (field_type)
640  {
642  {
643  const auto entity_it = _feature_maps[var_index].find(entity_id);
644 
645  if (entity_it != _feature_maps[var_index].end())
646  return entity_it->second; // + _region_offsets[var_index];
647  else
648  return -1;
649  }
650 
652  {
653  const auto entity_it = _var_index_maps[var_index].find(entity_id);
654 
655  if (entity_it != _var_index_maps[var_index].end())
656  return entity_it->second;
657  else
658  return -1;
659  }
660 
662  {
663  const auto entity_it = _ghosted_entity_ids.find(entity_id);
664 
665  if (entity_it != _ghosted_entity_ids.end())
666  return entity_it->second;
667  else
668  return -1;
669  }
670 
671  case FieldType::HALOS:
672  {
673  if (!use_default)
674  {
675  const auto entity_it = _halo_ids[var_index].find(entity_id);
676  if (entity_it != _halo_ids[var_index].end())
677  return entity_it->second;
678  }
679  else
680  {
681  // Showing halos in reverse order for backwards compatibility
682  for (auto map_num = _maps_size;
683  map_num-- /* don't compare greater than zero for unsigned */;)
684  {
685  const auto entity_it = _halo_ids[map_num].find(entity_id);
686 
687  if (entity_it != _halo_ids[map_num].end())
688  return entity_it->second;
689  }
690  }
691  return -1;
692  }
693 
694  case FieldType::CENTROID:
695  {
696  if (_periodic_node_map.size())
697  mooseDoOnce(mooseWarning(
698  "Centroids are not correct when using periodic boundaries, contact the MOOSE team"));
699 
700  // If this element contains the centroid of one of features, return one
701  const auto * elem_ptr = _mesh.elemPtr(entity_id);
702 
703  for (const auto & feature : _feature_sets)
704  {
705  if (feature._status == Status::INACTIVE)
706  continue;
707 
708  if (elem_ptr->contains_point(feature._centroid))
709  return 1;
710  }
711 
712  return 0;
713  }
714 
715  default:
716  return 0;
717  }
718 }
std::multimap< dof_id_type, dof_id_type > _periodic_node_map
The data structure which is a list of nodes that are constrained to other nodes based on the imposed ...
static const std::size_t invalid_size_t
std::map< dof_id_type, int > _ghosted_entity_ids
The map for holding reconstructed ghosted element information.
std::vector< std::map< dof_id_type, int > > _halo_ids
The data structure for looking up halos around features.
std::vector< FeatureData > _feature_sets
The data structure used to hold the globally unique features.
const std::size_t _maps_size
Convenience variable holding the size of all the datastructures size by the number of maps...
std::vector< std::map< dof_id_type, int > > _feature_maps
The feature maps contain the raw flooded node information and eventually the unique grain numbers...
MooseMesh & _mesh
A reference to the mesh.
std::vector< std::map< dof_id_type, int > > _var_index_maps
This map keeps track of which variables own which nodes.
const std::vector<FeatureData>& FeatureFloodCount::getFeatures ( ) const
inline

Return a constant reference to the vector of all discovered features.

Definition at line 290 of file FeatureFloodCount.h.

Referenced by GrainTracker::prepopulateState().

290 { return _feature_sets; }
std::vector< FeatureData > _feature_sets
The data structure used to hold the globally unique features.
unsigned int FeatureFloodCount::getFeatureVar ( unsigned int  feature_id) const
virtual

Returns the variable representing the passed in feature.

Reimplemented in GrainTracker, and FauxGrainTracker.

Definition at line 587 of file FeatureFloodCount.C.

Referenced by FeatureVolumeVectorPostprocessor::execute(), and GrainTracker::getFeatureVar().

588 {
589  // Some processors don't contain the largest feature id, in that case we just return invalid_id
590  if (feature_id >= _feature_id_to_local_index.size())
591  return invalid_id;
592 
593  auto local_index = _feature_id_to_local_index[feature_id];
594  if (local_index != invalid_size_t)
595  {
596  mooseAssert(local_index < _feature_sets.size(), "local_index out of bounds");
597  return _feature_sets[local_index]._status != Status::INACTIVE
598  ? _feature_sets[feature_id]._var_index
599  : invalid_id;
600  }
601 
602  return invalid_id;
603 }
static const std::size_t invalid_size_t
std::vector< FeatureData > _feature_sets
The data structure used to hold the globally unique features.
static const unsigned int invalid_id
std::vector< std::size_t > _feature_id_to_local_index
The vector recording the grain_id to local index (several indices will contain invalid_size_t) ...
Real FeatureFloodCount::getThreshold ( std::size_t  current_index) const
protectedvirtual

Return the starting comparison threshold to use when inspecting an entity during the flood stage.

Reimplemented in GrainTracker.

Definition at line 1072 of file FeatureFloodCount.C.

Referenced by getFeatures(), and isNewFeatureOrConnectedRegion().

1073 {
1074  return _step_threshold;
1075 }
std::size_t FeatureFloodCount::getTotalFeatureCount ( ) const
virtual

Returns the total feature count (active and inactive ids, useful for sizing vectors)

Since the FeatureFloodCount object doesn't maintain any information about features between invocations. The maximum id in use is simply the number of features.

Reimplemented in FauxGrainTracker, and GrainTracker.

Definition at line 576 of file FeatureFloodCount.C.

Referenced by FeatureVolumeVectorPostprocessor::execute(), and AverageGrainVolume::initialize().

577 {
583  return _feature_count;
584 }
unsigned int _feature_count
The number of features seen by this object (same as summing _feature_counts_per_map) ...
Real FeatureFloodCount::getValue ( )
overridevirtual

Reimplemented in FauxGrainTracker.

Definition at line 570 of file FeatureFloodCount.C.

571 {
572  return static_cast<Real>(_feature_count);
573 }
unsigned int _feature_count
The number of features seen by this object (same as summing _feature_counts_per_map) ...
const std::vector< unsigned int > & FeatureFloodCount::getVarToFeatureVector ( dof_id_type  elem_id) const
virtual

Returns a list of active unique feature ids for a particular element.

The vector is indexed by variable number with each entry containing either an invalid size_t type (no feature active at that location) or a feature id if the variable is non-zero at that location.

Reimplemented in GrainTracker, and FauxGrainTracker.

Definition at line 486 of file FeatureFloodCount.C.

Referenced by AverageGrainVolume::execute(), FeatureVolumeVectorPostprocessor::execute(), GrainTracker::getVarToFeatureVector(), and FeatureFloodCountAux::precalculateValue().

487 {
488  mooseDoOnce(if (!_compute_var_to_feature_map) mooseError(
489  "Please set \"compute_var_to_feature_map = true\" to use this interface method"));
490 
491  const auto pos = _entity_var_to_features.find(elem_id);
492  if (pos != _entity_var_to_features.end())
493  {
494  mooseAssert(pos->second.size() == _n_vars, "Variable to feature vector not sized properly");
495  return pos->second;
496  }
497  else
498  return _empty_var_to_features;
499 }
const std::size_t _n_vars
std::map< dof_id_type, std::vector< unsigned int > > _entity_var_to_features
std::vector< unsigned int > _empty_var_to_features
const bool _compute_var_to_feature_map
Indicates whether or not the var to feature map is populated.
void FeatureFloodCount::initialize ( )
overridevirtual

Reimplemented in PolycrystalUserObjectBase, FauxGrainTracker, and GrainTracker.

Definition at line 220 of file FeatureFloodCount.C.

Referenced by GrainTracker::initialize(), and PolycrystalUserObjectBase::initialize().

221 {
222  // Clear the feature marking maps and region counters and other data structures
223  for (auto map_num = decltype(_maps_size)(0); map_num < _maps_size; ++map_num)
224  {
225  _feature_maps[map_num].clear();
226  _partial_feature_sets[map_num].clear();
227 
228  if (_var_index_mode)
229  _var_index_maps[map_num].clear();
230 
231  _halo_ids[map_num].clear();
232  }
233 
234  _feature_sets.clear();
235 
236  // Calculate the thresholds for this iteration
239 
240  _ghosted_entity_ids.clear();
241 
242  // Reset the feature count and max local size
243  _feature_count = 0;
244 
245  _entity_var_to_features.clear();
246 
247  for (auto & map_ref : _entities_visited)
248  map_ref.clear();
249 }
const PostprocessorValue & _element_average_value
Average value of the domain which can optionally be used to find features in a field.
std::vector< std::set< dof_id_type > > _entities_visited
This variable keeps track of which nodes have been visited during execution.
std::map< dof_id_type, std::vector< unsigned int > > _entity_var_to_features
std::map< dof_id_type, int > _ghosted_entity_ids
The map for holding reconstructed ghosted element information.
std::vector< std::map< dof_id_type, int > > _halo_ids
The data structure for looking up halos around features.
const Real _connecting_threshold
The threshold above (or below) which neighboring entities are flooded (where regions can be extended ...
std::vector< FeatureData > _feature_sets
The data structure used to hold the globally unique features.
const std::size_t _maps_size
Convenience variable holding the size of all the datastructures size by the number of maps...
std::vector< std::map< dof_id_type, int > > _feature_maps
The feature maps contain the raw flooded node information and eventually the unique grain numbers...
unsigned int _feature_count
The number of features seen by this object (same as summing _feature_counts_per_map) ...
const Real _threshold
The threshold above (or below) where an entity may begin a new region (feature)
std::vector< std::list< FeatureData > > _partial_feature_sets
The data structure used to hold partial and communicated feature data.
const bool _var_index_mode
This variable is used to indicate whether the maps will contain unique region information or just the...
std::vector< std::map< dof_id_type, int > > _var_index_maps
This map keeps track of which variables own which nodes.
void FeatureFloodCount::initialSetup ( )
overridevirtual

Size the empty var to features vector to the number of coupled variables. This empty vector (but properly sized) vector is returned for elements that are queried but are not in the structure (which also shouldn't happen). The user is warned in this case but this helps avoid extra bounds checking in user code and avoids segfaults.

Reimplemented in PolycrystalUserObjectBase.

Definition at line 199 of file FeatureFloodCount.C.

Referenced by PolycrystalUserObjectBase::initialSetup().

200 {
201  // We need one map per coupled variable for normal runs to support overlapping features
202  _entities_visited.resize(_vars.size());
203 
204  // Get a pointer to the PeriodicBoundaries buried in libMesh
205  _pbs = _fe_problem.getNonlinearSystemBase().dofMap().get_periodic_boundaries();
206 
207  meshChanged();
208 
217 }
const std::size_t _n_vars
std::vector< std::set< dof_id_type > > _entities_visited
This variable keeps track of which nodes have been visited during execution.
std::vector< MooseVariable * > _vars
The vector of coupled in variables.
static const unsigned int invalid_id
PeriodicBoundaries * _pbs
A pointer to the periodic boundary constraints object.
std::vector< unsigned int > _empty_var_to_features
virtual void meshChanged() override
bool FeatureFloodCount::isElemental ( ) const
inline

Definition at line 102 of file FeatureFloodCount.h.

Referenced by FeatureFloodCountAux::FeatureFloodCountAux().

102 { return _is_elemental; }
bool _is_elemental
Determines if the flood counter is elements or not (nodes)
bool FeatureFloodCount::isNewFeatureOrConnectedRegion ( const DofObject *  dof_object,
std::size_t &  current_index,
FeatureData *&  feature,
Status status,
unsigned int &  new_id 
)
protectedvirtual

Method called during the recursive flood routine that should return whether or not the current entity is part of the current feature (if one is being explored), or if it's the start of a new feature.

If the value is only above the connecting threshold, it's still part of a feature but possibly part of one that we'll discard if there is never any starting threshold encountered.

Reimplemented in PolycrystalUserObjectBase.

Definition at line 1090 of file FeatureFloodCount.C.

Referenced by flood(), and getFeatures().

1095 {
1096  // Get the value of the current variable for the current entity
1097  Real entity_value;
1098  if (_is_elemental)
1099  {
1100  const Elem * elem = static_cast<const Elem *>(dof_object);
1101  std::vector<Point> centroid(1, elem->centroid());
1102  _subproblem.reinitElemPhys(elem, centroid, 0);
1103  entity_value = _vars[current_index]->sln()[0];
1104  }
1105  else
1106  entity_value = _vars[current_index]->getNodalValue(*static_cast<const Node *>(dof_object));
1107 
1108  // If the value compares against our starting threshold, this is definitely part of a feature
1109  // we'll keep
1110  if (compareValueWithThreshold(entity_value, getThreshold(current_index)))
1111  {
1112  Status * status_ptr = &status;
1113 
1114  if (feature)
1115  status_ptr = &feature->_status;
1116 
1117  // Update an existing feature's status or clear the flag on the passed in status
1118  *status_ptr &= ~Status::INACTIVE;
1119  return true;
1120  }
1121 
1126  return compareValueWithThreshold(entity_value, getConnectingThreshold(current_index));
1127 }
Status
This enumeration is used to indicate status of the grains in the _unique_grains data structure...
virtual Real getThreshold(std::size_t current_index) const
Return the starting comparison threshold to use when inspecting an entity during the flood stage...
std::vector< MooseVariable * > _vars
The vector of coupled in variables.
bool _is_elemental
Determines if the flood counter is elements or not (nodes)
bool compareValueWithThreshold(Real entity_value, Real threshold) const
This method is used to determine whether the current entity value is part of a feature or not...
virtual Real getConnectingThreshold(std::size_t current_index) const
Return the "connecting" comparison threshold to use when inspecting an entity during the flood stage...
void FeatureFloodCount::mergeSets ( )
protected

This routine is called on the master rank only and stitches together the partial feature pieces seen on any processor.

Insert the new entity at the end of the list so that it may be checked against all other partial features again.

Now remove both halves the merged features: it2 contains the "moved" feature cell just inserted at the back of the list, it1 contains the mostly empty other half. We have to be careful about the order in which these two elements are deleted. We delete it2 first since we don't care where its iterator points after the deletion. We are going to break out of this loop anyway. If we delete it1 first, it may end up pointing at the same location as it2 which after the second deletion would cause both of the iterators to be invalidated.

Now that the merges are complete we need to adjust the centroid, and halos. Additionally, To make several of the sorting and tracking algorithms more straightforward, we will move the features into a flat vector. Finally we can count the final number of features and find the max local index seen on any processor Note: This is all occurring on rank 0 only!

IMPORTANT: FeatureFloodCount::_feature_count is set on rank 0 at this point but we can't broadcast it here because this routine is not collective.

Definition at line 828 of file FeatureFloodCount.C.

Referenced by communicateAndMerge(), and getFeatures().

829 {
830  Moose::perf_log.push("mergeSets()", "FeatureFloodCount");
831 
832  // Since we gathered only on the root process, we only need to merge sets on the root process.
833  mooseAssert(_is_master, "mergeSets() should only be called on the root process");
834 
835  // Local variable used for sizing structures, it will be >= the actual number of features
836  for (auto map_num = decltype(_maps_size)(0); map_num < _maps_size; ++map_num)
837  {
838  for (auto it1 = _partial_feature_sets[map_num].begin();
839  it1 != _partial_feature_sets[map_num].end();
840  /* No increment on it1 */)
841  {
842  bool merge_occured = false;
843  for (auto it2 = _partial_feature_sets[map_num].begin();
844  it2 != _partial_feature_sets[map_num].end();
845  ++it2)
846  {
847  if (it1 != it2 && areFeaturesMergeable(*it1, *it2))
848  {
849  it2->merge(std::move(*it1));
850 
855  _partial_feature_sets[map_num].emplace_back(std::move(*it2));
856 
866  _partial_feature_sets[map_num].erase(it2);
867  it1 = _partial_feature_sets[map_num].erase(it1); // it1 is incremented here!
868 
869  // A merge occurred, this is used to determine whether or not we increment the outer
870  // iterator
871  merge_occured = true;
872 
873  // We need to start the list comparison over for the new it1 so break here
874  break;
875  }
876  } // it2 loop
877 
878  if (!merge_occured) // No merges so we need to manually increment the outer iterator
879  ++it1;
880 
881  } // it1 loop
882  } // map loop
883 
891  // Offset where the current set of features with the same variable id starts in the flat vector
892  unsigned int feature_offset = 0;
893  // Set the member feature count to zero and start counting the actual features
894  _feature_count = 0;
895 
896  for (auto map_num = decltype(_maps_size)(0); map_num < _maps_size; ++map_num)
897  {
898  std::set<dof_id_type> set_difference;
899  for (auto & feature : _partial_feature_sets[map_num])
900  {
901  // If after merging we still have an inactive feature, discard it
902  if (feature._status == Status::CLEAR)
903  {
904  // First we need to calculate the centroid now that we are doing merging all partial
905  // features
906  if (feature._vol_count != 0)
907  feature._centroid /= feature._vol_count;
908 
909  _feature_sets.emplace_back(std::move(feature));
910  ++_feature_count;
911  }
912  }
913 
914  // Record the feature numbers just for the current map
915  _feature_counts_per_map[map_num] = _feature_count - feature_offset;
916 
917  // Now update the running feature count so we can calculate the next map's contribution
918  feature_offset = _feature_count;
919 
920  // Clean up the "moved" objects
921  _partial_feature_sets[map_num].clear();
922  }
923 
929  Moose::perf_log.pop("mergeSets()", "FeatureFloodCount");
930 }
virtual bool areFeaturesMergeable(const FeatureData &f1, const FeatureData &f2) const
Method for determining whether two features are mergeable.
bool _is_master
Convenience variable for testing master rank.
std::vector< FeatureData > _feature_sets
The data structure used to hold the globally unique features.
const std::size_t _maps_size
Convenience variable holding the size of all the datastructures size by the number of maps...
unsigned int _feature_count
The number of features seen by this object (same as summing _feature_counts_per_map) ...
std::vector< std::list< FeatureData > > _partial_feature_sets
The data structure used to hold partial and communicated feature data.
std::vector< unsigned int > _feature_counts_per_map
The number of features seen by this object per map.
void FeatureFloodCount::meshChanged ( )
overridevirtual

We need to build a set containing all of the boundary entities to compare against. This will be elements for elemental flooding. Volumes for nodal flooding is not supported

Definition at line 257 of file FeatureFloodCount.C.

Referenced by initialSetup().

258 {
259  _point_locator = _mesh.getMesh().sub_point_locator();
260 
261  _mesh.buildPeriodicNodeMap(_periodic_node_map, _var_number, _pbs);
262 
263  // Build a new node to element map
264  _nodes_to_elem_map.clear();
265  MeshTools::build_nodes_to_elem_map(_mesh.getMesh(), _nodes_to_elem_map);
266 
272  _all_boundary_entity_ids.clear();
273  if (_is_elemental)
274  for (auto elem_it = _mesh.bndElemsBegin(), elem_end = _mesh.bndElemsEnd(); elem_it != elem_end;
275  ++elem_it)
276  _all_boundary_entity_ids.insert((*elem_it)->_elem->id());
277 }
std::multimap< dof_id_type, dof_id_type > _periodic_node_map
The data structure which is a list of nodes that are constrained to other nodes based on the imposed ...
bool _is_elemental
Determines if the flood counter is elements or not (nodes)
std::set< dof_id_type > _all_boundary_entity_ids
The set of entities on the boundary of the domain used for determining if features intersect any boun...
unsigned long _var_number
This variable is used to build the periodic node map.
std::vector< std::vector< const Elem * > > _nodes_to_elem_map
The data structure used to find neighboring elements give a node ID.
std::unique_ptr< PointLocatorBase > _point_locator
PeriodicBoundaries * _pbs
A pointer to the periodic boundary constraints object.
MooseMesh & _mesh
A reference to the mesh.
std::size_t FeatureFloodCount::numCoupledVars ( ) const
inline

Returns the number of coupled varaibles.

Definition at line 77 of file FeatureFloodCount.h.

77 { return _n_vars; }
const std::size_t _n_vars
void FeatureFloodCount::prepareDataForTransfer ( )
protected

This routine uses the local flooded data to build up the local feature data structures (_feature_sets).

This routine does not perform any communication so the _feature_sets data structure will only contain information from the local processor after calling this routine. Any existing data in the _feature_sets structure is destroyed by calling this routine.

_feature_sets layout: The outer vector is sized to one when _single_map_mode == true, otherwise it is sized for the number of coupled variables. The inner list represents the flooded regions (local only after this call but fully populated after parallel communication and stitching).

We need to adjust the halo markings before sending. We need to discard all of the local cell information but not any of the stitch region information. To do that we subtract off the ghosted cells from the local cells and use that in the set difference operation with the halo_ids.

Save off the min entity id present in the feature to uniquely identify the feature regardless of n_procs

Definition at line 721 of file FeatureFloodCount.C.

Referenced by communicateAndMerge(), and getFeatures().

722 {
723  MeshBase & mesh = _mesh.getMesh();
724 
725  std::set<dof_id_type> local_ids_no_ghost, set_difference;
726 
727  for (auto & list_ref : _partial_feature_sets)
728  {
729  for (auto & feature : list_ref)
730  {
731  // Now extend the bounding box by the halo region
732  if (_is_elemental)
733  feature.updateBBoxExtremes(mesh);
734  else
735  {
736  for (auto & halo_id : feature._halo_ids)
737  updateBBoxExtremesHelper(feature._bboxes[0], mesh.node(halo_id));
738  }
739 
746  std::set_difference(feature._local_ids.begin(),
747  feature._local_ids.end(),
748  feature._ghosted_ids.begin(),
749  feature._ghosted_ids.end(),
750  std::insert_iterator<std::set<dof_id_type>>(local_ids_no_ghost,
751  local_ids_no_ghost.begin()));
752 
753  std::set_difference(
754  feature._halo_ids.begin(),
755  feature._halo_ids.end(),
756  local_ids_no_ghost.begin(),
757  local_ids_no_ghost.end(),
758  std::insert_iterator<std::set<dof_id_type>>(set_difference, set_difference.begin()));
759  feature._halo_ids.swap(set_difference);
760  local_ids_no_ghost.clear();
761  set_difference.clear();
762 
763  mooseAssert(!feature._local_ids.empty(), "local entity ids cannot be empty");
764 
769  feature._min_entity_id = *feature._local_ids.begin();
770 
771  // Periodic node ids
773  }
774  }
775 }
void appendPeriodicNeighborNodes(FeatureData &feature) const
This routine adds the periodic node information to our data structure prior to packing the data this ...
bool _is_elemental
Determines if the flood counter is elements or not (nodes)
void updateBBoxExtremesHelper(MeshTools::BoundingBox &bbox, const Point &node)
std::vector< std::list< FeatureData > > _partial_feature_sets
The data structure used to hold partial and communicated feature data.
MooseMesh & _mesh
A reference to the mesh.
void FeatureFloodCount::scatterAndUpdateRanks ( )
protected

Calls buildLocalToGlobalIndices to build the individual local to global indicies for each rank and scatters that information to all ranks.

Finally, the non-master ranks update their own data structures to reflect the global mappings.

On non-root processors we can't maintain the full _feature_sets data structure since we don't have all of the global information. We'll move the items from the partial feature sets into a flat structure maintaining order and update the internal IDs with the proper global ID.

Important: Make sure we clear the local status if we received a valid global index for this feature. It's possible that we have a status of INVALID on the local processor because there was never any starting threshold found. However, the root processor wouldn't have sent an index if it didn't find a starting threshold connected to our local piece.

Definition at line 502 of file FeatureFloodCount.C.

Referenced by GrainTracker::assignGrains(), finalize(), getFeatures(), and GrainTracker::trackGrains().

503 {
504  // local to global map (one per processor)
505  std::vector<int> counts;
506  std::vector<std::size_t> local_to_global_all;
507  if (_is_master)
508  buildLocalToGlobalIndices(local_to_global_all, counts);
509 
510  // Scatter local_to_global indices to all processors and store in class member variable
511  _communicator.scatter(local_to_global_all, counts, _local_to_global_feature_map);
512 
513  std::size_t largest_global_index = std::numeric_limits<std::size_t>::lowest();
514  if (!_is_master)
515  {
517 
524  for (auto & list_ref : _partial_feature_sets)
525  {
526  for (auto & feature : list_ref)
527  {
528  mooseAssert(feature._orig_ids.size() == 1, "feature._orig_ids length doesn't make sense");
529 
530  auto global_index = FeatureFloodCount::invalid_size_t;
531  auto local_index = feature._orig_ids.begin()->second;
532 
533  if (local_index < _local_to_global_feature_map.size())
534  global_index = _local_to_global_feature_map[local_index];
535 
536  if (global_index != FeatureFloodCount::invalid_size_t)
537  {
538  if (global_index > largest_global_index)
539  largest_global_index = global_index;
540 
541  // Set the correct global index
542  feature._id = global_index;
543 
551  feature._status &= ~Status::INACTIVE;
552 
553  // Move the feature into the correct place
554  _feature_sets[local_index] = std::move(feature);
555  }
556  }
557  }
558  }
559  else
560  {
561  for (auto global_index : local_to_global_all)
562  if (global_index != FeatureFloodCount::invalid_size_t && global_index > largest_global_index)
563  largest_global_index = global_index;
564  }
565 
566  buildFeatureIdToLocalIndices(largest_global_index);
567 }
static const std::size_t invalid_size_t
Status
This enumeration is used to indicate status of the grains in the _unique_grains data structure...
bool _is_master
Convenience variable for testing master rank.
std::vector< FeatureData > _feature_sets
The data structure used to hold the globally unique features.
virtual void buildLocalToGlobalIndices(std::vector< std::size_t > &local_to_global_all, std::vector< int > &counts) const
This routine populates a stacked vector of local to global indices per rank and the associated count ...
std::vector< std::list< FeatureData > > _partial_feature_sets
The data structure used to hold partial and communicated feature data.
void buildFeatureIdToLocalIndices(unsigned int max_id)
This method builds a lookup map for retrieving the right local feature (by index) given a global inde...
std::vector< std::size_t > _local_to_global_feature_map
The vector recording the local to global feature indices.
void FeatureFloodCount::serialize ( std::string &  serialized_buffer)
protected

These routines packs/unpack the _feature_map data into a structure suitable for parallel communication operations.

See the comments in these routines for the exact data structure layout.

Call the MOOSE serialization routines to serialize this processor's data. Note: The _partial_feature_sets data structure will be empty for all other processors

Definition at line 778 of file FeatureFloodCount.C.

Referenced by communicateAndMerge(), and getFeatures().

779 {
780  // stream for serializing the _partial_feature_sets data structure to a byte stream
781  std::ostringstream oss;
782 
787  dataStore(oss, _partial_feature_sets, this);
788 
789  // Populate the passed in string pointer with the string stream's buffer contents
790  serialized_buffer.assign(oss.str());
791 }
void dataStore(std::ostream &stream, FeatureFloodCount::FeatureData &feature, void *context)
std::vector< std::list< FeatureData > > _partial_feature_sets
The data structure used to hold partial and communicated feature data.
template<class InputIterator >
static bool FeatureFloodCount::setsIntersect ( InputIterator  first1,
InputIterator  last1,
InputIterator  first2,
InputIterator  last2 
)
inlinestaticprotected

This method detects whether two sets intersect without building a result set.

It exits as soon as any intersection is detected.

Definition at line 482 of file FeatureFloodCount.h.

Referenced by FeatureFloodCount::FeatureData::ghostedIntersect(), FeatureFloodCount::FeatureData::halosIntersect(), and FeatureFloodCount::FeatureData::periodicBoundariesIntersect().

486  {
487  while (first1 != last1 && first2 != last2)
488  {
489  if (*first1 == *first2)
490  return true;
491 
492  if (*first1 < *first2)
493  ++first1;
494  else if (*first1 > *first2)
495  ++first2;
496  }
497  return false;
498  }
void FeatureFloodCount::sortAndLabel ( )
protected

Sort and assign ids to features based on their position in the container after sorting.

Perform a sort to give a parallel unique sorting to the identified features. We use the "min_entity_id" inside each feature to assign it's position in the sorted vector.

Sanity check. Now that we've sorted the flattened vector of features we need to make sure that the counts vector still lines up appropriately with each feature's _var_index.

Definition at line 361 of file FeatureFloodCount.C.

Referenced by GrainTracker::assignGrains(), finalize(), and getFeatures().

362 {
363  mooseAssert(_is_master, "sortAndLabel can only be called on the master");
364 
370  std::sort(_feature_sets.begin(), _feature_sets.end());
371 
372 #ifndef NDEBUG
373 
378  unsigned int feature_offset = 0;
379  for (auto map_num = beginIndex(_feature_counts_per_map); map_num < _maps_size; ++map_num)
380  {
381  // Skip empty map checks
382  if (_feature_counts_per_map[map_num] == 0)
383  continue;
384 
385  // Check the begin and end of the current range
386  auto range_front = feature_offset;
387  auto range_back = feature_offset + _feature_counts_per_map[map_num] - 1;
388 
389  mooseAssert(range_front <= range_back && range_back < _feature_count,
390  "Indexing error in feature sets");
391 
392  if (!_single_map_mode && (_feature_sets[range_front]._var_index != map_num ||
393  _feature_sets[range_back]._var_index != map_num))
394  mooseError("Error in _feature_sets sorting, map index: ", map_num);
395 
396  feature_offset += _feature_counts_per_map[map_num];
397  }
398 #endif
399 
400  // Label the features with an ID based on the sorting (processor number independent value)
401  for (auto i = beginIndex(_feature_sets); i < _feature_sets.size(); ++i)
402  if (_feature_sets[i]._id == invalid_id)
403  _feature_sets[i]._id = i;
404 }
bool _is_master
Convenience variable for testing master rank.
std::vector< FeatureData > _feature_sets
The data structure used to hold the globally unique features.
const std::size_t _maps_size
Convenience variable holding the size of all the datastructures size by the number of maps...
static const unsigned int invalid_id
const bool _single_map_mode
This variable is used to indicate whether or not multiple maps are used during flooding.
unsigned int _feature_count
The number of features seen by this object (same as summing _feature_counts_per_map) ...
std::vector< unsigned int > _feature_counts_per_map
The number of features seen by this object per map.
void FeatureFloodCount::updateFieldInfo ( )
protectedvirtual

This method is used to populate any of the data structures used for storing field data (nodal or elemental).

It is called at the end of finalize and can make use of any of the data structures created during the execution of this postprocessor.

Reimplemented in GrainTracker.

Definition at line 939 of file FeatureFloodCount.C.

Referenced by finalize(), and getFeatures().

940 {
941  for (auto i = beginIndex(_feature_sets); i < _feature_sets.size(); ++i)
942  {
943  auto & feature = _feature_sets[i];
944 
945  // If the developer has requested _condense_map_info we'll make sure we only update the zeroth
946  // map
947  auto map_index = (_single_map_mode || _condense_map_info) ? decltype(feature._var_index)(0)
948  : feature._var_index;
949 
950  // Loop over the entity ids of this feature and update our local map
951  for (auto entity : feature._local_ids)
952  {
953  _feature_maps[map_index][entity] = static_cast<int>(feature._id);
954 
955  if (_var_index_mode)
956  _var_index_maps[map_index][entity] = feature._var_index;
957 
958  // Fill in the data structure that keeps track of all features per elem
960  {
961  auto insert_pair = moose_try_emplace(
962  _entity_var_to_features, entity, std::vector<unsigned int>(_n_vars, invalid_id));
963  auto & vec_ref = insert_pair.first->second;
964  vec_ref[feature._var_index] = feature._id;
965  }
966  }
967 
968  if (_compute_halo_maps)
969  // Loop over the halo ids to update cells with halo information
970  for (auto entity : feature._halo_ids)
971  _halo_ids[map_index][entity] = static_cast<int>(feature._id);
972 
973  // Loop over the ghosted ids to update cells with ghost information
974  for (auto entity : feature._ghosted_ids)
975  _ghosted_entity_ids[entity] = 1;
976 
977  // TODO: Fixme
978  if (!_global_numbering)
979  mooseError("Local numbering currently disabled");
980  }
981 }
const std::size_t _n_vars
const bool _condense_map_info
std::map< dof_id_type, std::vector< unsigned int > > _entity_var_to_features
std::map< dof_id_type, int > _ghosted_entity_ids
The map for holding reconstructed ghosted element information.
std::vector< std::map< dof_id_type, int > > _halo_ids
The data structure for looking up halos around features.
std::vector< FeatureData > _feature_sets
The data structure used to hold the globally unique features.
const bool _global_numbering
This variable is used to indicate whether or not we identify features with unique numbers on multiple...
static const unsigned int invalid_id
std::vector< std::map< dof_id_type, int > > _feature_maps
The feature maps contain the raw flooded node information and eventually the unique grain numbers...
const bool _single_map_mode
This variable is used to indicate whether or not multiple maps are used during flooding.
const bool _compute_halo_maps
Indicates whether or not to communicate halo map information with all ranks.
const bool _compute_var_to_feature_map
Indicates whether or not the var to feature map is populated.
const bool _var_index_mode
This variable is used to indicate whether the maps will contain unique region information or just the...
std::vector< std::map< dof_id_type, int > > _var_index_maps
This map keeps track of which variables own which nodes.
void FeatureFloodCount::updateRegionOffsets ( )
protected

This routine updates the _region_offsets variable which is useful for quickly determining the proper global number for a feature when using multimap mode.

Referenced by getFeatures().

void FeatureFloodCount::visitElementalNeighbors ( const Elem *  elem,
std::size_t  current_index,
FeatureData feature,
bool  expand_halos_only,
bool  disjoint_only 
)
protected

Retrieve only the active neighbors for each side of this element, append them to the list of active neighbors

If the current element (passed into this method) doesn't have a connected neighbor but does have a topological neighbor, this might be a new disjoint region that we'll need to represent with a separate bounding box. To find out for sure, we'll need see if the new neighbors are present in any of the halo or disjoint halo sets. If they are not present, this is a new region.

Definition at line 1246 of file FeatureFloodCount.C.

Referenced by expandEdgeHalos(), flood(), and getFeatures().

1251 {
1252  mooseAssert(elem, "Elem is NULL");
1253 
1254  std::vector<const Elem *> all_active_neighbors;
1255  MeshBase & mesh = _mesh.getMesh();
1256 
1257  // Loop over all neighbors (at the the same level as the current element)
1258  // Loop over all neighbors (at the the same level as the current element)
1259  for (auto i = decltype(elem->n_neighbors())(0); i < elem->n_neighbors(); ++i)
1260  {
1261  const Elem * neighbor_ancestor = nullptr;
1262  bool topological_neighbor = false;
1263 
1268  neighbor_ancestor = elem->neighbor(i);
1269  if (neighbor_ancestor)
1270  neighbor_ancestor->active_family_tree_by_neighbor(all_active_neighbors, elem, false);
1271  else // if (expand_halos_only /*&& feature->_periodic_nodes.empty()*/)
1272  {
1273  neighbor_ancestor = elem->topological_neighbor(i, mesh, *_point_locator, _pbs);
1274 
1282  if (neighbor_ancestor)
1283  {
1284  neighbor_ancestor->active_family_tree_by_topological_neighbor(
1285  all_active_neighbors, elem, mesh, *_point_locator, _pbs, false);
1286 
1287  topological_neighbor = true;
1288 
1289  // for (const auto neighbor : all_active_neighbors)
1290  // feature->_disjoint_halo_ids.insert(neighbor->id());
1291  }
1292  }
1293 
1294  visitNeighborsHelper(elem,
1295  all_active_neighbors,
1296  current_index,
1297  feature,
1298  expand_halos_only,
1299  topological_neighbor,
1300  disjoint_only);
1301 
1302  all_active_neighbors.clear();
1303  }
1304 }
void visitNeighborsHelper(const T *curr_entity, std::vector< const T * > neighbor_entities, std::size_t current_index, FeatureData *feature, bool expand_halos_only, bool topological_neighbor, bool disjoint_only)
The actual logic for visiting neighbors is abstracted out here.
std::unique_ptr< PointLocatorBase > _point_locator
PeriodicBoundaries * _pbs
A pointer to the periodic boundary constraints object.
MooseMesh & _mesh
A reference to the mesh.
template<typename T >
void FeatureFloodCount::visitNeighborsHelper ( const T *  curr_entity,
std::vector< const T * >  neighbor_entities,
std::size_t  current_index,
FeatureData feature,
bool  expand_halos_only,
bool  topological_neighbor,
bool  disjoint_only 
)
protected

The actual logic for visiting neighbors is abstracted out here.

This method is templated to handle the Nodal and Elemental cases together.

Only recurse where we own this entity and it's a topologically connected entity. We shouldn't even attempt to flood to the periodic boundary because we won't have solution information and if we are using DistributedMesh we probably won't have geometric information either.

When we only recurse on entities we own, we can never get more than one away from a local entity which should be in the ghosted zone.

Premark neighboring entities with a halo mark. These entities may or may not end up being part of the feature. We will not update the _entities_visited data structure here.

Definition at line 1323 of file FeatureFloodCount.C.

Referenced by getFeatures(), visitElementalNeighbors(), and visitNodalNeighbors().

1330 {
1331  // Loop over all active element neighbors
1332  for (const auto neighbor : neighbor_entities)
1333  {
1334  if (neighbor)
1335  {
1336  if (expand_halos_only)
1337  {
1338  if (topological_neighbor || disjoint_only)
1339  feature->_disjoint_halo_ids.insert(neighbor->id());
1340  else
1341  feature->_halo_ids.insert(neighbor->id());
1342  }
1343  else
1344  {
1345  auto my_processor_id = processor_id();
1346 
1347  if (!topological_neighbor && neighbor->processor_id() != my_processor_id)
1348  feature->_ghosted_ids.insert(curr_entity->id());
1349 
1359  if (curr_entity->processor_id() == my_processor_id)
1360  {
1367  if (topological_neighbor || disjoint_only)
1368  feature->_disjoint_halo_ids.insert(neighbor->id());
1369  else
1370  {
1371  feature->_halo_ids.insert(neighbor->id());
1372 
1373  flood(neighbor, current_index, feature);
1374  }
1375  }
1376  }
1377  }
1378  }
1379 }
bool flood(const DofObject *dof_object, std::size_t current_index, FeatureData *feature)
This method will "mark" all entities on neighboring elements that are above the supplied threshold...
void FeatureFloodCount::visitNodalNeighbors ( const Node *  node,
std::size_t  current_index,
FeatureData feature,
bool  expand_halos_only 
)
protected

These two routines are utility routines used by the flood routine and by derived classes for visiting neighbors.

Since the logic is different for the elemental versus nodal case it's easier to split them up.

Definition at line 1307 of file FeatureFloodCount.C.

Referenced by expandEdgeHalos(), flood(), and getFeatures().

1311 {
1312  mooseAssert(node, "Node is NULL");
1313 
1314  std::vector<const Node *> all_active_neighbors;
1315  MeshTools::find_nodal_neighbors(_mesh.getMesh(), *node, _nodes_to_elem_map, all_active_neighbors);
1316 
1318  node, all_active_neighbors, current_index, feature, expand_halos_only, false, false);
1319 }
void visitNeighborsHelper(const T *curr_entity, std::vector< const T * > neighbor_entities, std::size_t current_index, FeatureData *feature, bool expand_halos_only, bool topological_neighbor, bool disjoint_only)
The actual logic for visiting neighbors is abstracted out here.
std::vector< std::vector< const Elem * > > _nodes_to_elem_map
The data structure used to find neighboring elements give a node ID.
MooseMesh & _mesh
A reference to the mesh.

Member Data Documentation

std::set<dof_id_type> FeatureFloodCount::_all_boundary_entity_ids
protected

The set of entities on the boundary of the domain used for determining if features intersect any boundary.

Definition at line 637 of file FeatureFloodCount.h.

Referenced by flood(), and meshChanged().

const bool FeatureFloodCount::_compute_halo_maps
protected

Indicates whether or not to communicate halo map information with all ranks.

Definition at line 540 of file FeatureFloodCount.h.

Referenced by GrainTracker::communicateHaloMap(), GrainTracker::updateFieldInfo(), and updateFieldInfo().

const bool FeatureFloodCount::_compute_var_to_feature_map
protected

Indicates whether or not the var to feature map is populated.

Definition at line 543 of file FeatureFloodCount.h.

Referenced by getVarToFeatureVector(), GrainTracker::updateFieldInfo(), and updateFieldInfo().

const bool FeatureFloodCount::_condense_map_info
protected

Definition at line 529 of file FeatureFloodCount.h.

Referenced by GrainTracker::updateFieldInfo(), and updateFieldInfo().

const Real FeatureFloodCount::_connecting_threshold
protected

The threshold above (or below) which neighboring entities are flooded (where regions can be extended but not started)

Definition at line 513 of file FeatureFloodCount.h.

Referenced by initialize().

const PostprocessorValue& FeatureFloodCount::_element_average_value
protected

Average value of the domain which can optionally be used to find features in a field.

Definition at line 618 of file FeatureFloodCount.h.

Referenced by initialize().

std::vector<unsigned int> FeatureFloodCount::_empty_var_to_features
protected

Definition at line 641 of file FeatureFloodCount.h.

Referenced by getVarToFeatureVector(), and initialSetup().

std::vector<std::set<dof_id_type> > FeatureFloodCount::_entities_visited
protected

This variable keeps track of which nodes have been visited during execution.

We don't use the _feature_map for this since we don't want to explicitly store data for all the unmarked nodes in a serialized datastructures. This keeps our overhead down since this variable never needs to be communicated.

Definition at line 567 of file FeatureFloodCount.h.

Referenced by PolycrystalUserObjectBase::execute(), flood(), initialize(), initialSetup(), and PolycrystalUserObjectBase::isNewFeatureOrConnectedRegion().

std::map<dof_id_type, std::vector<unsigned int> > FeatureFloodCount::_entity_var_to_features
protected
unsigned int FeatureFloodCount::_feature_count
protected
std::vector<unsigned int> FeatureFloodCount::_feature_counts_per_map
protected

The number of features seen by this object per map.

Definition at line 581 of file FeatureFloodCount.h.

Referenced by mergeSets(), sortAndLabel(), and GrainTracker::trackGrains().

std::vector<std::size_t> FeatureFloodCount::_feature_id_to_local_index
protected
std::vector<std::map<dof_id_type, int> > FeatureFloodCount::_feature_maps
protected

The feature maps contain the raw flooded node information and eventually the unique grain numbers.

We have a vector of them so we can create one per variable if that level of detail is desired.

Definition at line 604 of file FeatureFloodCount.h.

Referenced by getEntityValue(), initialize(), GrainTracker::updateFieldInfo(), and updateFieldInfo().

std::vector<FeatureData> FeatureFloodCount::_feature_sets
protected
std::map<dof_id_type, int> FeatureFloodCount::_ghosted_entity_ids
protected

The map for holding reconstructed ghosted element information.

Definition at line 621 of file FeatureFloodCount.h.

Referenced by getEntityValue(), initialize(), GrainTracker::updateFieldInfo(), and updateFieldInfo().

const bool FeatureFloodCount::_global_numbering
protected

This variable is used to indicate whether or not we identify features with unique numbers on multiple maps.

Definition at line 533 of file FeatureFloodCount.h.

Referenced by updateFieldInfo().

std::vector<std::map<dof_id_type, int> > FeatureFloodCount::_halo_ids
protected

The data structure for looking up halos around features.

The outer vector is for splitting out the information per variable. The inner map holds the actual halo information

Definition at line 627 of file FeatureFloodCount.h.

Referenced by FeatureFloodCount::FeatureData::clear(), GrainTracker::communicateHaloMap(), getEntityValue(), FeatureFloodCount::FeatureData::halosIntersect(), initialize(), FeatureFloodCount::FeatureData::merge(), FeatureFloodCount::FeatureData::updateBBoxExtremes(), GrainTracker::updateFieldInfo(), and updateFieldInfo().

bool FeatureFloodCount::_is_elemental
protected
bool FeatureFloodCount::_is_master
protected
std::vector<std::size_t> FeatureFloodCount::_local_to_global_feature_map
protected

The vector recording the local to global feature indices.

Definition at line 607 of file FeatureFloodCount.h.

Referenced by scatterAndUpdateRanks().

const std::size_t FeatureFloodCount::_maps_size
protected

Convenience variable holding the size of all the datastructures size by the number of maps.

Definition at line 556 of file FeatureFloodCount.h.

Referenced by FeatureFloodCount(), getEntityValue(), initialize(), mergeSets(), sortAndLabel(), GrainTracker::trackGrains(), and GrainTracker::updateFieldInfo().

MooseMesh& FeatureFloodCount::_mesh
protected
const processor_id_type FeatureFloodCount::_n_procs
protected

Convenience variable holding the number of processors in this simulation.

Definition at line 559 of file FeatureFloodCount.h.

Referenced by buildLocalToGlobalIndices(), and GrainTracker::communicateHaloMap().

const std::size_t FeatureFloodCount::_n_vars
protected
std::vector<std::vector<const Elem *> > FeatureFloodCount::_nodes_to_elem_map
protected

The data structure used to find neighboring elements give a node ID.

Definition at line 578 of file FeatureFloodCount.h.

Referenced by meshChanged(), and visitNodalNeighbors().

std::vector<std::list<FeatureData> > FeatureFloodCount::_partial_feature_sets
protected

The data structure used to hold partial and communicated feature data.

The data structure mirrors that found in _feature_sets, but contains one additional vector indexed by processor id

Definition at line 591 of file FeatureFloodCount.h.

Referenced by deserialize(), expandEdgeHalos(), expandPointHalos(), flood(), initialize(), mergeSets(), prepareDataForTransfer(), GrainTracker::prepopulateState(), scatterAndUpdateRanks(), and serialize().

PeriodicBoundaries* FeatureFloodCount::_pbs
protected

A pointer to the periodic boundary constraints object.

Definition at line 613 of file FeatureFloodCount.h.

Referenced by initialSetup(), meshChanged(), and visitElementalNeighbors().

std::multimap<dof_id_type, dof_id_type> FeatureFloodCount::_periodic_node_map
protected

The data structure which is a list of nodes that are constrained to other nodes based on the imposed periodic boundary conditions.

Definition at line 633 of file FeatureFloodCount.h.

Referenced by appendPeriodicNeighborNodes(), FauxGrainTracker::getEntityValue(), getEntityValue(), and meshChanged().

std::unique_ptr<PointLocatorBase> FeatureFloodCount::_point_locator
protected

Definition at line 615 of file FeatureFloodCount.h.

Referenced by meshChanged(), and visitElementalNeighbors().

const bool FeatureFloodCount::_single_map_mode
protected

This variable is used to indicate whether or not multiple maps are used during flooding.

Definition at line 527 of file FeatureFloodCount.h.

Referenced by flood(), PolycrystalUserObjectBase::PolycrystalUserObjectBase(), sortAndLabel(), GrainTracker::updateFieldInfo(), and updateFieldInfo().

Real FeatureFloodCount::_step_connecting_threshold
protected

Definition at line 514 of file FeatureFloodCount.h.

Referenced by getConnectingThreshold(), and initialize().

Real FeatureFloodCount::_step_threshold
protected

Definition at line 509 of file FeatureFloodCount.h.

Referenced by GrainTracker::getThreshold(), getThreshold(), and initialize().

const Real FeatureFloodCount::_threshold
protected

The threshold above (or below) where an entity may begin a new region (feature)

Definition at line 508 of file FeatureFloodCount.h.

Referenced by FauxGrainTracker::execute(), and initialize().

const bool FeatureFloodCount::_use_less_than_threshold_comparison
protected

Use less-than when comparing values against the threshold value.

True by default. If false, then greater-than comparison is used instead.

Definition at line 550 of file FeatureFloodCount.h.

Referenced by compareValueWithThreshold(), and FauxGrainTracker::execute().

std::vector<std::map<dof_id_type, int> > FeatureFloodCount::_var_index_maps
protected

This map keeps track of which variables own which nodes.

We need a vector of them for multimap mode where multiple variables can own a single mode.

Note: This map is only populated when "show_var_coloring" is set to true.

Definition at line 575 of file FeatureFloodCount.h.

Referenced by FeatureFloodCount(), getEntityValue(), initialize(), GrainTracker::updateFieldInfo(), and updateFieldInfo().

const bool FeatureFloodCount::_var_index_mode
protected

This variable is used to indicate whether the maps will contain unique region information or just the variable numbers owning those regions.

Definition at line 537 of file FeatureFloodCount.h.

Referenced by FeatureFloodCount(), initialize(), GrainTracker::updateFieldInfo(), and updateFieldInfo().

unsigned long FeatureFloodCount::_var_number
protected

This variable is used to build the periodic node map.

Assumption: We are going to assume that either all variables are periodic or none are. This assumption can be relaxed at a later time if necessary.

Definition at line 524 of file FeatureFloodCount.h.

Referenced by GrainTracker::centroidRegionDistance(), and meshChanged().

std::vector<MooseVariable *> FeatureFloodCount::_vars
protected
const unsigned int FeatureFloodCount::invalid_id = std::numeric_limits<unsigned int>::max()
static
const std::size_t FeatureFloodCount::invalid_size_t = std::numeric_limits<std::size_t>::max()
static

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