www.mooseframework.org
GrainTracker.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
3 /* */
4 /* All contents are licensed under LGPL V2.1 */
5 /* See LICENSE for full restrictions */
6 /****************************************************************/
7 
8 #include "GrainTracker.h"
9 
10 // MOOSE includes
12 #include "GeneratedMesh.h"
13 #include "MooseMesh.h"
14 #include "MooseVariable.h"
15 #include "NonlinearSystem.h"
16 
17 #include "libmesh/periodic_boundary_base.h"
18 
19 // C++ includes
20 #include <algorithm>
21 #include <limits>
22 #include <numeric>
23 
24 template <>
25 void
26 dataStore(std::ostream & stream, GrainTracker::PartialFeatureData & feature, void * context)
27 {
28  storeHelper(stream, feature.intersects_boundary, context);
29  storeHelper(stream, feature.id, context);
30  storeHelper(stream, feature.centroid, context);
31  storeHelper(stream, feature.status, context);
32 }
33 
34 template <>
35 void
36 dataLoad(std::istream & stream, GrainTracker::PartialFeatureData & feature, void * context)
37 {
38  loadHelper(stream, feature.intersects_boundary, context);
39  loadHelper(stream, feature.id, context);
40  loadHelper(stream, feature.centroid, context);
41  loadHelper(stream, feature.status, context);
42 }
43 
44 template <>
45 InputParameters
47 {
48  InputParameters params = validParams<FeatureFloodCount>();
50  params.addClassDescription("Grain Tracker object for running reduced order parameter simulations "
51  "without grain coalescence.");
52 
53  return params;
54 }
55 
56 GrainTracker::GrainTracker(const InputParameters & parameters)
57  : FeatureFloodCount(parameters),
59  _tracking_step(getParam<int>("tracking_step")),
60  _halo_level(getParam<unsigned int>("halo_level")),
61  _n_reserve_ops(getParam<unsigned short>("reserve_op")),
62  _reserve_op_index(_n_reserve_ops <= _n_vars ? _n_vars - _n_reserve_ops : 0),
63  _reserve_op_threshold(getParam<Real>("reserve_op_threshold")),
64  _remap(getParam<bool>("remap_grains")),
65  _nl(_fe_problem.getNonlinearSystemBase()),
66  _feature_sets_old(declareRestartableData<std::vector<FeatureData>>("unique_grains")),
67  _poly_ic_uo(parameters.isParamValid("polycrystal_ic_uo")
68  ? &getUserObject<PolycrystalUserObjectBase>("polycrystal_ic_uo")
69  : nullptr),
70  _first_time(true),
71  _error_on_grain_creation(getParam<bool>("error_on_grain_creation")),
72  _reserve_grain_first_index(0),
73  _old_max_grain_id(0),
74  _max_curr_grain_id(0),
75  _is_transient(_subproblem.isTransient())
76 {
77  if (_tracking_step > 0 && _poly_ic_uo)
78  mooseError("Can't start tracking after the initial condition when using a polycrystal_ic_uo");
79 }
80 
82 
83 Real
84 GrainTracker::getEntityValue(dof_id_type entity_id,
85  FieldType field_type,
86  std::size_t var_index) const
87 {
88  if (_t_step < _tracking_step)
89  return 0;
90 
91  return FeatureFloodCount::getEntityValue(entity_id, field_type, var_index);
92 }
93 
94 const std::vector<unsigned int> &
95 GrainTracker::getVarToFeatureVector(dof_id_type elem_id) const
96 {
98 }
99 
100 unsigned int
101 GrainTracker::getFeatureVar(unsigned int feature_id) const
102 {
103  return FeatureFloodCount::getFeatureVar(feature_id);
104 }
105 
106 std::size_t
108 {
109  // Note: This value is parallel consistent, see FeatureFloodCount::communicateAndMerge()
110  return _feature_count;
111 }
112 
113 std::size_t
115 {
116  // Note: This value is parallel consistent, see assignGrains()/trackGrains()
117  return _max_curr_grain_id + 1;
118 }
119 
120 Point
121 GrainTracker::getGrainCentroid(unsigned int grain_id) const
122 {
123  mooseAssert(grain_id < _feature_id_to_local_index.size(), "Grain ID out of bounds");
124  auto grain_index = _feature_id_to_local_index[grain_id];
125 
126  if (grain_index != invalid_size_t)
127  {
128  mooseAssert(_feature_id_to_local_index[grain_id] < _feature_sets.size(),
129  "Grain index out of bounds");
130  // Note: This value is parallel consistent, see GrainTracker::broadcastAndUpdateGrainData()
131  return _feature_sets[_feature_id_to_local_index[grain_id]]._centroid;
132  }
133 
134  // Inactive grain
135  return Point();
136 }
137 
138 bool
139 GrainTracker::doesFeatureIntersectBoundary(unsigned int feature_id) const
140 {
141  // TODO: This data structure may need to be turned into a Multimap
142  mooseAssert(feature_id < _feature_id_to_local_index.size(), "Grain ID out of bounds");
143 
144  auto feature_index = _feature_id_to_local_index[feature_id];
145  if (feature_index != invalid_size_t)
146  {
147  mooseAssert(feature_index < _feature_sets.size(), "Grain index out of bounds");
148  return _feature_sets[feature_index]._intersects_boundary;
149  }
150 
151  return false;
152 }
153 
154 void
156 {
157  // Don't track grains if the current simulation step is before the specified tracking step
158  if (_t_step < _tracking_step)
159  return;
160 
166  if (!_first_time)
168 
170 }
171 
172 void
174 {
175  // Don't track grains if the current simulation step is before the specified tracking step
176  if (_t_step < _tracking_step)
177  return;
178 
179  if (_poly_ic_uo && _first_time)
180  return;
181 
182  Moose::perf_log.push("execute()", "GrainTracker");
184  Moose::perf_log.pop("execute()", "GrainTracker");
185 }
186 
187 Real
188 GrainTracker::getThreshold(std::size_t var_index) const
189 {
190  // If we are inspecting a reserve op parameter, we need to make sure
191  // that there is an entity above the reserve_op threshold before
192  // starting the flood of the feature.
193  if (var_index >= _reserve_op_index)
194  return _reserve_op_threshold;
195  else
196  return _step_threshold;
197 }
198 
199 void
201 {
202  mooseAssert(_first_time, "This method should only be called on the first invocation");
203 
204  _feature_sets.clear();
205 
211  if (_is_master)
212  {
213  const auto & features = ffc_object.getFeatures();
214  for (auto & feature : features)
215  _feature_sets.emplace_back(feature.duplicate());
216 
217  _feature_count = _feature_sets.size();
218  }
219  else
220  {
221  const auto & features = ffc_object.getFeatures();
222  _partial_feature_sets[0].clear();
223  for (auto & feature : features)
224  _partial_feature_sets[0].emplace_back(feature.duplicate());
225  }
226 
227  // Make sure that feature count is communicated to all ranks
228  _communicator.broadcast(_feature_count);
229 }
230 
231 void
233 {
239  // Don't track grains if the current simulation step is before the specified tracking step
240  if (_t_step < _tracking_step)
241  return;
242 
243  Moose::perf_log.push("finalize()", "GrainTracker");
244 
245  // Expand the depth of the halos around all grains
246  auto num_halo_layers = _halo_level >= 1
247  ? _halo_level - 1
248  : 0; // The first level of halos already exists so subtract one
249 
250  if (_poly_ic_uo && _first_time)
252  else
253  {
254  expandEdgeHalos(num_halo_layers);
255 
256  // Build up the grain map on the root processor
258  }
259 
263  Moose::perf_log.push("trackGrains()", "GrainTracker");
264  if (_first_time)
265  assignGrains();
266  else
267  trackGrains();
268  Moose::perf_log.pop("trackGrains()", "GrainTracker");
269  _console << "Finished inside of trackGrains" << std::endl;
270 
275 
279  Moose::perf_log.push("remapGrains()", "GrainTracker");
280  if (_remap)
281  remapGrains();
282  Moose::perf_log.pop("remapGrains()", "GrainTracker");
283 
284  updateFieldInfo();
285  _console << "Finished inside of updateFieldInfo" << std::endl;
286 
287  // Set the first time flag false here (after all methods of finalize() have completed)
288  _first_time = false;
289 
290  // TODO: Release non essential memory
291 
292  _console << "Finished inside of GrainTracker" << std::endl;
293  Moose::perf_log.pop("finalize()", "GrainTracker");
294 }
295 
296 void
298 {
299  std::vector<PartialFeatureData> root_feature_data;
300  std::vector<std::string> send_buffer(1), recv_buffer;
301 
302  if (_is_master)
303  {
304  root_feature_data.reserve(_feature_sets.size());
305 
306  // Populate a subset of the information in a small data structure
307  std::transform(_feature_sets.begin(),
308  _feature_sets.end(),
309  std::back_inserter(root_feature_data),
310  [](FeatureData & feature) {
311  PartialFeatureData partial_feature;
312  partial_feature.intersects_boundary = feature._intersects_boundary;
313  partial_feature.id = feature._id;
314  partial_feature.centroid = feature._centroid;
315  partial_feature.status = feature._status;
316  return partial_feature;
317  });
318 
319  std::ostringstream oss;
320  dataStore(oss, root_feature_data, this);
321  send_buffer[0].assign(oss.str());
322  }
323 
324  // Broadcast the data to all ranks
325  _communicator.broadcast_packed_range((void *)(nullptr),
326  send_buffer.begin(),
327  send_buffer.end(),
328  (void *)(nullptr),
329  std::back_inserter(recv_buffer));
330 
331  // Unpack and update
332  if (!_is_master)
333  {
334  std::istringstream iss;
335  iss.str(recv_buffer[0]);
336  iss.clear();
337 
338  dataLoad(iss, root_feature_data, this);
339 
340  for (const auto & partial_data : root_feature_data)
341  {
342  // See if this processor has a record of this grain
343  if (partial_data.id < _feature_id_to_local_index.size() &&
344  _feature_id_to_local_index[partial_data.id] != invalid_size_t)
345  {
346  auto & grain = _feature_sets[_feature_id_to_local_index[partial_data.id]];
347  grain._intersects_boundary = partial_data.intersects_boundary;
348  grain._centroid = partial_data.centroid;
349  if (partial_data.status == Status::INACTIVE)
350  grain._status = Status::INACTIVE;
351  }
352  }
353  }
354 }
355 
356 void
358 {
359  mooseAssert(_first_time, "assignGrains may only be called on the first tracking step");
360 
366  if (_is_master)
367  {
368  mooseAssert(!_feature_sets.empty(), "Feature sets empty!");
369 
370  // Find the largest grain ID, this requires sorting if the ID is not already set
371  sortAndLabel();
373 
374  for (auto & grain : _feature_sets)
375  grain._status = Status::MARKED; // Mark the grain
376 
377  // Set up the first reserve grain index based on the largest grain ID
379  } // is_master
380 
381  /*************************************************************
382  ****************** COLLECTIVE WORK SECTION ******************
383  *************************************************************/
384 
385  // Make IDs on all non-master ranks consistent
387 
388  // Build up an id to index map
389  _communicator.broadcast(_max_curr_grain_id);
391 
392  // Now trigger the newGrainCreated() callback on all ranks
393  for (auto new_id = decltype(_max_curr_grain_id)(0); new_id <= _max_curr_grain_id; ++new_id)
394  newGrainCreated(new_id);
395 }
396 
397 void
399 {
400  mooseAssert(!_first_time, "Track grains may only be called when _tracking_step > _t_step");
401 
402  // Used to track indices for which to trigger the new grain callback on (used on all ranks)
404 
409  if (_is_master)
410  {
411  // Reset Status on active unique grains
412  std::vector<unsigned int> map_sizes(_maps_size);
413  for (auto & grain : _feature_sets_old)
414  {
415  if (grain._status != Status::INACTIVE)
416  {
417  grain._status = Status::CLEAR;
418  map_sizes[grain._var_index]++;
419  }
420  }
421 
422  // Print out stats on overall tracking changes per var_index
423  for (auto map_num = decltype(_maps_size)(0); map_num < _maps_size; ++map_num)
424  {
425  _console << "\nGrains active index " << map_num << ": " << map_sizes[map_num] << " -> "
426  << _feature_counts_per_map[map_num];
427  if (map_sizes[map_num] > _feature_counts_per_map[map_num])
428  _console << "--";
429  else if (map_sizes[map_num] < _feature_counts_per_map[map_num])
430  _console << "++";
431  }
432  _console << '\n' << std::endl;
433 
440  std::vector<std::size_t> new_grain_index_to_existing_grain_index(_feature_sets.size(),
442 
443  for (auto old_grain_index = beginIndex(_feature_sets_old);
444  old_grain_index < _feature_sets_old.size();
445  ++old_grain_index)
446  {
447  auto & old_grain = _feature_sets_old[old_grain_index];
448 
449  if (old_grain._status == Status::INACTIVE) // Don't try to find matches for inactive grains
450  continue;
451 
452  std::size_t closest_match_index = invalid_size_t;
453  Real min_centroid_diff = std::numeric_limits<Real>::max();
454 
460  // clang-format off
461  auto start_it =
462  std::lower_bound(_feature_sets.begin(), _feature_sets.end(), old_grain._var_index,
463  [](const FeatureData & item, std::size_t var_index)
464  {
465  return item._var_index < var_index;
466  });
467  // clang-format on
468 
469  // We only need to examine grains that have matching variable indices
470  for (decltype(_feature_sets.size()) new_grain_index =
471  std::distance(_feature_sets.begin(), start_it);
472  new_grain_index < _feature_sets.size() &&
473  _feature_sets[new_grain_index]._var_index == old_grain._var_index;
474  ++new_grain_index)
475  {
476  auto & new_grain = _feature_sets[new_grain_index];
477 
482  if (new_grain.boundingBoxesIntersect(old_grain))
483  {
484  Real curr_centroid_diff = centroidRegionDistance(old_grain._bboxes, new_grain._bboxes);
485  if (curr_centroid_diff <= min_centroid_diff)
486  {
487  closest_match_index = new_grain_index;
488  min_centroid_diff = curr_centroid_diff;
489  }
490  }
491  }
492 
493  // found a match
494  if (closest_match_index != invalid_size_t)
495  {
502  auto curr_index = new_grain_index_to_existing_grain_index[closest_match_index];
503  if (curr_index != invalid_size_t)
504  {
505  // The new feature being competed for
506  auto & new_grain = _feature_sets[closest_match_index];
507 
508  // The other old grain competing to match up to the same new grain
509  auto & other_old_grain = _feature_sets_old[curr_index];
510 
511  auto centroid_diff1 = centroidRegionDistance(new_grain._bboxes, old_grain._bboxes);
512  auto centroid_diff2 = centroidRegionDistance(new_grain._bboxes, other_old_grain._bboxes);
513 
514  auto & inactive_grain = (centroid_diff1 < centroid_diff2) ? other_old_grain : old_grain;
515 
516  inactive_grain._status = Status::INACTIVE;
517  _console << COLOR_GREEN << "Marking Grain " << inactive_grain._id
518  << " as INACTIVE (variable index: " << inactive_grain._var_index << ")\n"
519  << COLOR_DEFAULT << inactive_grain;
520 
526  if (&inactive_grain == &other_old_grain)
527  new_grain_index_to_existing_grain_index[closest_match_index] = old_grain_index;
528  }
529  else
530  new_grain_index_to_existing_grain_index[closest_match_index] = old_grain_index;
531  }
532  }
533 
534  // Mark all resolved grain matches
535  for (auto new_index = beginIndex(new_grain_index_to_existing_grain_index);
536  new_index < new_grain_index_to_existing_grain_index.size();
537  ++new_index)
538  {
539  auto curr_index = new_grain_index_to_existing_grain_index[new_index];
540 
541  // This may be a new grain, we'll handle that case below
542  if (curr_index == invalid_size_t)
543  continue;
544 
545  mooseAssert(_feature_sets_old[curr_index]._id != invalid_id,
546  "Invalid ID in old grain structure");
547 
548  _feature_sets[new_index]._id = _feature_sets_old[curr_index]._id; // Transfer ID
549  _feature_sets[new_index]._status = Status::MARKED; // Mark the status in the new set
550  _feature_sets_old[curr_index]._status = Status::MARKED; // Mark the status in the old set
551  }
552 
567  // Case 1 (new grains in _feature_sets):
568  for (auto grain_num = beginIndex(_feature_sets); grain_num < _feature_sets.size(); ++grain_num)
569  {
570  auto & grain = _feature_sets[grain_num];
571 
572  // New Grain
573  if (grain._status == Status::CLEAR)
574  {
593  // clang-format off
594  auto start_it =
595  std::lower_bound(_feature_sets.begin(), _feature_sets.end(), grain._var_index,
596  [](const FeatureData & item, std::size_t var_index)
597  {
598  return item._var_index < var_index;
599  });
600  // clang-format on
601 
602  // Loop over matching variable indices
603  for (decltype(_feature_sets.size()) new_grain_index =
604  std::distance(_feature_sets.begin(), start_it);
605  new_grain_index < _feature_sets.size() &&
606  _feature_sets[new_grain_index]._var_index == grain._var_index;
607  ++new_grain_index)
608  {
609  auto & other_grain = _feature_sets[new_grain_index];
610 
611  // Splitting grain?
612  if (grain_num != new_grain_index && // Make sure indices aren't pointing at the same grain
613  other_grain._status == Status::MARKED && // and that the other grain is indeed marked
614  other_grain.boundingBoxesIntersect(grain) && // and the bboxes intersect
615  other_grain.halosIntersect(grain)) // and the halos also intersect
616  // TODO: Inspect combined volume and see if it's "close" to the expected value
617  {
618  grain._id = other_grain._id; // Set the duplicate ID
619  grain._status = Status::MARKED; // Mark it
620  _console << COLOR_YELLOW << "Split Grain Detected "
621  << " (variable index: " << grain._var_index << ")\n"
622  << COLOR_DEFAULT << grain << other_grain;
623  }
624  }
625 
626  // Must be a nucleating grain (status is still not set)
627  if (grain._status == Status::CLEAR)
628  {
629  auto new_index = getNextUniqueID();
630  grain._id = new_index; // Set the ID
631  grain._status = Status::MARKED; // Mark it
632  }
633  }
634  }
635 
636  // Case 2 (inactive grains in _feature_sets_old)
637  for (auto & grain : _feature_sets_old)
638  {
639  if (grain._status == Status::CLEAR)
640  {
641  grain._status = Status::INACTIVE;
642  _console << COLOR_GREEN << "Marking Grain " << grain._id
643  << " as INACTIVE (variable index: " << grain._var_index << ")\n"
644  << COLOR_DEFAULT << grain;
645  }
646  }
647  } // is_master
648 
649  /*************************************************************
650  ****************** COLLECTIVE WORK SECTION ******************
651  *************************************************************/
652 
653  // Make IDs on all non-master ranks consistent
655 
656  // Build up an id to index map
657  _communicator.broadcast(_max_curr_grain_id);
659 
664  {
665  for (auto new_id = _old_max_grain_id + 1; new_id <= _max_curr_grain_id; ++new_id)
666  {
667  // Don't trigger the callback on the reserve IDs
669  {
670  // See if we've been instructed to terminate with an error
672  mooseError(
673  "Error: New grain detected and \"error_on_new_grain_creation\" is set to true");
674  else
675  newGrainCreated(new_id);
676  }
677  }
678  }
679 }
680 
681 void
682 GrainTracker::newGrainCreated(unsigned int new_grain_id)
683 {
684  if (!_first_time && _is_master)
685  {
686  mooseAssert(new_grain_id < _feature_id_to_local_index.size(), "new_grain_id is out of bounds");
687  auto grain_index = _feature_id_to_local_index[new_grain_id];
688  mooseAssert(grain_index != invalid_size_t && grain_index < _feature_sets.size(),
689  "new_grain_id appears to be invalid");
690 
691  const auto & grain = _feature_sets[grain_index];
692  _console << COLOR_YELLOW
693  << "\n*****************************************************************************"
694  << "\nCouldn't find a matching grain while working on variable index: "
695  << grain._var_index << "\nCreating new unique grain: " << new_grain_id << '\n'
696  << grain
697  << "\n*****************************************************************************\n"
698  << COLOR_DEFAULT;
699  }
700 }
701 
702 std::vector<unsigned int>
704 {
705  std::vector<unsigned int> new_ids(_max_curr_grain_id - _old_max_grain_id);
706  auto new_id = _old_max_grain_id + 1;
707 
708  // Generate the new ids
709  std::iota(new_ids.begin(), new_ids.end(), new_id);
710 
711  return new_ids;
712 }
713 
714 void
716 {
717  // Don't remap grains if the current simulation step is before the specified tracking step
718  if (_t_step < _tracking_step)
719  return;
720 
721  _console << "Running remap Grains" << std::endl;
722 
729  std::map<unsigned int, std::size_t> grain_id_to_new_var;
730 
731  // Items are added to this list when split grains are found
732  std::list<std::pair<std::size_t, std::size_t>> split_pairs;
733 
742  if (_is_master)
743  {
744  // Build the map to detect difference in _var_index mappings after the remap operation
745  std::map<unsigned int, std::size_t> grain_id_to_existing_var_index;
746  for (auto & grain : _feature_sets)
747  {
748  // Unmark the grain so it can be used in the remap loop
749  grain._status = Status::CLEAR;
750 
751  grain_id_to_existing_var_index[grain._id] = grain._var_index;
752  }
753 
754  // Make sure that all split pieces of any grain are on the same OP
755  for (auto i = beginIndex(_feature_sets); i < _feature_sets.size(); ++i)
756  {
757  auto & grain1 = _feature_sets[i];
758 
759  for (auto j = beginIndex(_feature_sets); j < _feature_sets.size(); ++j)
760  {
761  auto & grain2 = _feature_sets[j];
762  if (i == j)
763  continue;
764 
765  // The first condition below is there to prevent symmetric checks (duplicate values)
766  if (i < j && grain1._id == grain2._id)
767  {
768  split_pairs.push_front(std::make_pair(i, j));
769  if (grain1._var_index != grain2._var_index)
770  {
771  _console << COLOR_YELLOW << "Split Grain (#" << grain1._id
772  << ") detected on unmatched OPs (" << grain1._var_index << ", "
773  << grain2._var_index << ") attempting to remap to " << grain1._var_index
774  << ".\n"
775  << COLOR_DEFAULT;
776 
781  grain1._var_index = grain2._var_index;
782  grain1._status |= Status::DIRTY;
783  }
784  }
785  }
786  }
787 
791  bool any_grains_remapped = false;
792  bool grains_remapped;
793  do
794  {
795  grains_remapped = false;
796  for (auto & grain1 : _feature_sets)
797  {
798  // We need to remap any grains represented on any variable index above the cuttoff
799  if (grain1._var_index >= _reserve_op_index)
800  {
801  _console << COLOR_YELLOW << "\nGrain #" << grain1._id
802  << " detected on a reserved order parameter #" << grain1._var_index
803  << ", remapping to another variable\n"
804  << COLOR_DEFAULT;
805 
806  for (auto max = decltype(_max_renumbering_recursion)(0);
808  ++max)
809  if (max < _max_renumbering_recursion)
810  {
811  if (attemptGrainRenumber(grain1, 0, max))
812  break;
813  }
814  else if (!attemptGrainRenumber(grain1, 0, max))
815  {
816  _console << std::flush;
817  mooseError(COLOR_RED,
818  "Unable to find any suitable order parameters for remapping."
819  " Perhaps you need more op variables?\n\n",
820  COLOR_DEFAULT);
821  }
822 
823  grains_remapped = true;
824  }
825 
826  for (auto & grain2 : _feature_sets)
827  {
828  // Don't compare a grain with itself and don't try to remap inactive grains
829  if (&grain1 == &grain2)
830  continue;
831 
832  if (grain1._var_index == grain2._var_index && // grains represented by same variable?
833  grain1._id != grain2._id && // are they part of different grains?
834  grain1.boundingBoxesIntersect(grain2) && // do bboxes intersect (coarse level)?
835  grain1.halosIntersect(grain2)) // do they actually overlap (fine level)?
836  {
837  _console << COLOR_YELLOW << "\nGrain #" << grain1._id << " intersects Grain #"
838  << grain2._id << " (variable index: " << grain1._var_index << ")\n"
839  << COLOR_DEFAULT;
840 
841  for (auto max = decltype(_max_renumbering_recursion)(0);
843  ++max)
844  if (max < _max_renumbering_recursion)
845  {
846  if (attemptGrainRenumber(grain1, 0, max))
847  break;
848  }
849  else if (!attemptGrainRenumber(grain1, 0, max) &&
850  !attemptGrainRenumber(grain2, 0, max))
851  {
852  _console << std::flush;
853  mooseError(COLOR_RED,
854  "Unable to find any suitable order parameters for remapping."
855  " Perhaps you need more op variables?\n\n",
856  COLOR_DEFAULT);
857  }
858 
859  grains_remapped = true;
860  }
861  }
862  }
863  any_grains_remapped |= grains_remapped;
864  } while (grains_remapped);
865 
866  // Verify that split grains are still intact
867  for (auto & split_pair : split_pairs)
868  if (_feature_sets[split_pair.first]._var_index != _feature_sets[split_pair.first]._var_index)
869  mooseError("Split grain remapped - This case is currently not handled");
870 
876  for (auto & grain : _feature_sets)
877  {
878  mooseAssert(grain_id_to_existing_var_index.find(grain._id) !=
879  grain_id_to_existing_var_index.end(),
880  "Missing unique ID");
881 
882  auto old_var_index = grain_id_to_existing_var_index[grain._id];
883 
884  if (old_var_index != grain._var_index)
885  {
886  mooseAssert(static_cast<bool>(grain._status & Status::DIRTY), "grain status is incorrect");
887 
888  grain_id_to_new_var.emplace_hint(
889  grain_id_to_new_var.end(),
890  std::pair<unsigned int, std::size_t>(grain._id, grain._var_index));
891 
900  grain._var_index = old_var_index;
901  // Clear the DIRTY status as well for consistency
902  grain._status &= ~Status::DIRTY;
903  }
904  }
905 
906  if (!grain_id_to_new_var.empty())
907  {
908  _console << "\nFinal remapping tally:\n";
909  for (const auto & remap_pair : grain_id_to_new_var)
910  _console << "Grain #" << remap_pair.first << " var_index "
911  << grain_id_to_existing_var_index[remap_pair.first] << " -> " << remap_pair.second
912  << '\n';
913  _console << "Communicating swaps with remaining processors..." << std::endl;
914  }
915  } // root processor
916 
917  // Communicate the std::map to all ranks
918  _communicator.broadcast(grain_id_to_new_var);
919 
920  // Perform swaps if any occurred
921  if (!grain_id_to_new_var.empty())
922  {
923  // Cache for holding values during swaps
924  std::vector<std::map<Node *, CacheValues>> cache(_n_vars);
925 
926  // Perform the actual swaps on all processors
927  for (auto & grain : _feature_sets)
928  {
929  // See if this grain was remapped
930  auto new_var_it = grain_id_to_new_var.find(grain._id);
931  if (new_var_it != grain_id_to_new_var.end())
932  swapSolutionValues(grain, new_var_it->second, cache, RemapCacheMode::FILL);
933  }
934 
935  for (auto & grain : _feature_sets)
936  {
937  // See if this grain was remapped
938  auto new_var_it = grain_id_to_new_var.find(grain._id);
939  if (new_var_it != grain_id_to_new_var.end())
940  swapSolutionValues(grain, new_var_it->second, cache, RemapCacheMode::USE);
941  }
942 
943  _nl.solution().close();
944  _nl.solutionOld().close();
945  _nl.solutionOlder().close();
946 
947  _fe_problem.getNonlinearSystemBase().system().update();
948 
949  _console << "Swaps complete" << std::endl;
950  }
951 }
952 
953 void
955  std::vector<std::list<GrainDistance>> & min_distances)
956 {
980  for (auto i = beginIndex(_feature_sets); i < _feature_sets.size(); ++i)
981  {
982  auto & other_grain = _feature_sets[i];
983 
984  if (other_grain._var_index == grain._var_index || other_grain._var_index >= _reserve_op_index)
985  continue;
986 
987  auto target_var_index = other_grain._var_index;
988  auto target_grain_index = i;
989  auto target_grain_id = other_grain._id;
990 
991  Real curr_bbox_diff = boundingRegionDistance(grain._bboxes, other_grain._bboxes);
992 
993  GrainDistance grain_distance_obj(
994  curr_bbox_diff, target_var_index, target_grain_index, target_grain_id);
995 
996  // To handle touching halos we penalize the top pick each time we see another
997  if (curr_bbox_diff == -1.0 && !min_distances[target_var_index].empty())
998  {
999  Real last_distance = min_distances[target_var_index].begin()->_distance;
1000  if (last_distance < 0)
1001  grain_distance_obj._distance += last_distance;
1002  }
1003 
1004  // Insertion sort into a list
1005  auto insert_it = min_distances[target_var_index].begin();
1006  while (insert_it != min_distances[target_var_index].end() && !(grain_distance_obj < *insert_it))
1007  ++insert_it;
1008  min_distances[target_var_index].insert(insert_it, grain_distance_obj);
1009  }
1010 
1016  for (auto var_index = beginIndex(_vars); var_index < _reserve_op_index; ++var_index)
1017  {
1018  // Don't put an entry in for matching variable indices (i.e. we can't remap to ourselves)
1019  if (grain._var_index == var_index)
1020  continue;
1021 
1022  if (min_distances[var_index].empty())
1023  min_distances[var_index].emplace_front(std::numeric_limits<Real>::max(), var_index);
1024  }
1025 }
1026 
1027 bool
1028 GrainTracker::attemptGrainRenumber(FeatureData & grain, unsigned int depth, unsigned int max_depth)
1029 {
1030  // End the recursion of our breadth first search
1031  if (depth > max_depth)
1032  return false;
1033 
1034  std::size_t curr_var_index = grain._var_index;
1035 
1036  std::vector<std::map<Node *, CacheValues>> cache;
1037 
1038  std::vector<std::list<GrainDistance>> min_distances(_vars.size());
1039 
1045  computeMinDistancesFromGrain(grain, min_distances);
1046 
1057  // clang-format off
1058  std::sort(min_distances.begin(), min_distances.end(),
1059  [](const std::list<GrainDistance> & lhs, const std::list<GrainDistance> & rhs)
1060  {
1061  // Sort lists in reverse order (largest distance first)
1062  // These empty cases are here to make this comparison stable
1063  if (lhs.empty())
1064  return false;
1065  else if (rhs.empty())
1066  return true;
1067  else
1068  return lhs.begin()->_distance > rhs.begin()->_distance;
1069  });
1070  // clang-format on
1071 
1072  for (auto & list_ref : min_distances)
1073  {
1074  const auto target_it = list_ref.begin();
1075  if (target_it == list_ref.end())
1076  continue;
1077 
1078  // If the distance is positive we can just remap and be done
1079  if (target_it->_distance > 0)
1080  {
1081  _console << COLOR_GREEN << "- Depth " << depth << ": Remapping grain #" << grain._id
1082  << " from variable index " << curr_var_index << " to " << target_it->_var_index;
1083  if (target_it->_distance == std::numeric_limits<Real>::max())
1084  _console << " which currently contains zero grains." << COLOR_DEFAULT;
1085  else
1086  _console << " whose closest grain (#" << target_it->_grain_id << ") is at a distance of "
1087  << target_it->_distance << "\n"
1088  << COLOR_DEFAULT;
1089 
1090  grain._status |= Status::DIRTY;
1091  grain._var_index = target_it->_var_index;
1092  return true;
1093  }
1094 
1095  // If the distance isn't positive we just need to make sure that none of the grains represented
1096  // by the target variable index would intersect this one if we were to remap
1097  auto next_target_it = target_it;
1098  bool intersection_hit = false;
1099  std::ostringstream oss;
1100  while (!intersection_hit && next_target_it != list_ref.end())
1101  {
1102  if (next_target_it->_distance > 0)
1103  break;
1104 
1105  mooseAssert(next_target_it->_grain_index < _feature_sets.size(),
1106  "Error in indexing target grain in attemptGrainRenumber");
1107  FeatureData & next_target_grain = _feature_sets[next_target_it->_grain_index];
1108 
1109  // If any grains touch we're done here
1110  if (grain.halosIntersect(next_target_grain))
1111  intersection_hit = true;
1112  else
1113  oss << " #" << next_target_it->_grain_id;
1114 
1115  ++next_target_it;
1116  }
1117 
1118  if (!intersection_hit)
1119  {
1120  _console << COLOR_GREEN << "- Depth " << depth << ": Remapping grain #" << grain._id
1121  << " from variable index " << curr_var_index << " to " << target_it->_var_index
1122  << " whose closest grain:" << oss.str()
1123  << " is inside our bounding box but whose halo(s) are not touching.\n"
1124  << COLOR_DEFAULT;
1125 
1126  grain._status |= Status::DIRTY;
1127  grain._var_index = target_it->_var_index;
1128  return true;
1129  }
1130 
1131  // If we reach this part of the loop, there is no simple renumbering that can be done.
1132  mooseAssert(target_it->_grain_index < _feature_sets.size(),
1133  "Error in indexing target grain in attemptGrainRenumber");
1134  FeatureData & target_grain = _feature_sets[target_it->_grain_index];
1135 
1143  if (target_it->_distance < -1)
1144  return false;
1145 
1146  // Make sure this grain isn't marked. If it is, we can't recurse here
1147  if ((target_grain._status & Status::MARKED) == Status::MARKED)
1148  return false;
1149 
1155  grain._var_index = target_it->_var_index;
1156  grain._status |= Status::MARKED;
1157  if (attemptGrainRenumber(target_grain, depth + 1, max_depth))
1158  {
1159  // SUCCESS!
1160  _console << COLOR_GREEN << "- Depth " << depth << ": Remapping grain #" << grain._id
1161  << " from variable index " << curr_var_index << " to " << target_it->_var_index
1162  << '\n'
1163  << COLOR_DEFAULT;
1164 
1165  // Now we need to mark the grain as DIRTY since the recursion succeeded.
1166  grain._status |= Status::DIRTY;
1167  return true;
1168  }
1169  else
1170  // FAILURE, We need to set our var index back after failed recursive step
1171  grain._var_index = curr_var_index;
1172 
1173  // ALWAYS "unmark" (or clear the MARKED status) after recursion so it can be used by other remap
1174  // operations
1175  grain._status &= ~Status::MARKED;
1176  }
1177 
1178  return false;
1179 }
1180 
1181 void
1183  std::size_t new_var_index,
1184  std::vector<std::map<Node *, CacheValues>> & cache,
1185  RemapCacheMode cache_mode)
1186 {
1187  MeshBase & mesh = _mesh.getMesh();
1188 
1189  // Remap the grain
1190  std::set<Node *> updated_nodes_tmp; // Used only in the elemental case
1191  for (auto entity : grain._local_ids)
1192  {
1193  if (_is_elemental)
1194  {
1195  Elem * elem = mesh.query_elem(entity);
1196  if (!elem)
1197  continue;
1198 
1199  for (unsigned int i = 0; i < elem->n_nodes(); ++i)
1200  {
1201  Node * curr_node = elem->get_node(i);
1202  if (updated_nodes_tmp.find(curr_node) == updated_nodes_tmp.end())
1203  {
1204  // cache this node so we don't attempt to remap it again within this loop
1205  updated_nodes_tmp.insert(curr_node);
1206  swapSolutionValuesHelper(curr_node, grain._var_index, new_var_index, cache, cache_mode);
1207  }
1208  }
1209  }
1210  else
1212  mesh.query_node_ptr(entity), grain._var_index, new_var_index, cache, cache_mode);
1213  }
1214 
1215  // Update the variable index in the unique grain datastructure after swaps are complete
1216  if (cache_mode == RemapCacheMode::USE || cache_mode == RemapCacheMode::BYPASS)
1217  grain._var_index = new_var_index;
1218 }
1219 
1220 void
1222  std::size_t curr_var_index,
1223  std::size_t new_var_index,
1224  std::vector<std::map<Node *, CacheValues>> & cache,
1225  RemapCacheMode cache_mode)
1226 {
1227  if (curr_node && curr_node->processor_id() == processor_id())
1228  {
1229  // Reinit the node so we can get and set values of the solution here
1230  _subproblem.reinitNode(curr_node, 0);
1231 
1232  // Local variables to hold values being transferred
1233  Real current, old = 0, older = 0;
1234  // Retrieve the value either from the old variable or cache
1235  if (cache_mode == RemapCacheMode::FILL || cache_mode == RemapCacheMode::BYPASS)
1236  {
1237  current = _vars[curr_var_index]->nodalSln()[0];
1238  if (_is_transient)
1239  {
1240  old = _vars[curr_var_index]->nodalSlnOld()[0];
1241  older = _vars[curr_var_index]->nodalSlnOlder()[0];
1242  }
1243  }
1244  else // USE
1245  {
1246  const auto cache_it = cache[curr_var_index].find(curr_node);
1247  mooseAssert(cache_it != cache[curr_var_index].end(), "Error in cache");
1248  current = cache_it->second.current;
1249  old = cache_it->second.old;
1250  older = cache_it->second.older;
1251  }
1252 
1253  // Cache the value or use it!
1254  if (cache_mode == RemapCacheMode::FILL)
1255  {
1256  cache[curr_var_index][curr_node].current = current;
1257  cache[curr_var_index][curr_node].old = old;
1258  cache[curr_var_index][curr_node].older = older;
1259  }
1260  else // USE or BYPASS
1261  {
1262  const auto & dof_index = _vars[new_var_index]->nodalDofIndex();
1263 
1264  // Transfer this solution from the old to the current
1265  _nl.solution().set(dof_index, current);
1266  if (_is_transient)
1267  {
1268  _nl.solutionOld().set(dof_index, old);
1269  _nl.solutionOlder().set(dof_index, older);
1270  }
1271  }
1272 
1285  if (cache_mode == RemapCacheMode::FILL || cache_mode == RemapCacheMode::BYPASS)
1286  {
1287  const auto & dof_index = _vars[curr_var_index]->nodalDofIndex();
1288 
1289  // Set the DOF for the current variable to zero
1290  _nl.solution().set(dof_index, 0.0);
1291  if (_is_transient)
1292  {
1293  _nl.solutionOld().set(dof_index, 0.0);
1294  _nl.solutionOlder().set(dof_index, 0.0);
1295  }
1296  }
1297  }
1298 }
1299 
1300 void
1302 {
1303  for (auto map_num = decltype(_maps_size)(0); map_num < _maps_size; ++map_num)
1304  _feature_maps[map_num].clear();
1305 
1306  std::map<dof_id_type, Real> tmp_map;
1307  MeshBase & mesh = _mesh.getMesh();
1308 
1309  for (const auto & grain : _feature_sets)
1310  {
1311  std::size_t curr_var = grain._var_index;
1312  std::size_t map_index = (_single_map_mode || _condense_map_info) ? 0 : curr_var;
1313 
1314  for (auto entity : grain._local_ids)
1315  {
1316  // Highest variable value at this entity wins
1317  Real entity_value = std::numeric_limits<Real>::lowest();
1318  if (_is_elemental)
1319  {
1320  const Elem * elem = mesh.elem(entity);
1321  std::vector<Point> centroid(1, elem->centroid());
1322  if (_poly_ic_uo && _first_time)
1323  {
1324  entity_value = _poly_ic_uo->getVariableValue(grain._var_index, centroid[0]);
1325  }
1326  else
1327  {
1328  _fe_problem.reinitElemPhys(elem, centroid, 0);
1329  entity_value = _vars[curr_var]->sln()[0];
1330  }
1331  }
1332  else
1333  {
1334  Node & node = mesh.node(entity);
1335  entity_value = _vars[curr_var]->getNodalValue(node);
1336  }
1337 
1338  if (entity_value != std::numeric_limits<Real>::lowest() &&
1339  (tmp_map.find(entity) == tmp_map.end() || entity_value > tmp_map[entity]))
1340  {
1341  mooseAssert(grain._id != invalid_id, "Missing Grain ID");
1342  _feature_maps[map_index][entity] = grain._id;
1343 
1344  if (_var_index_mode)
1345  _var_index_maps[map_index][entity] = grain._var_index;
1346 
1347  tmp_map[entity] = entity_value;
1348  }
1349 
1351  {
1352  auto insert_pair = moose_try_emplace(
1353  _entity_var_to_features, entity, std::vector<unsigned int>(_n_vars, invalid_id));
1354  auto & vec_ref = insert_pair.first->second;
1355 
1356  if (insert_pair.second)
1357  {
1358  // insert the reserve op numbers (if appropriate)
1359  for (auto reserve_index = decltype(_n_reserve_ops)(0); reserve_index < _n_reserve_ops;
1360  ++reserve_index)
1361  vec_ref[reserve_index] = _reserve_grain_first_index + reserve_index;
1362  }
1363  vec_ref[grain._var_index] = grain._id;
1364  }
1365  }
1366 
1367  if (_compute_halo_maps)
1368  for (auto entity : grain._halo_ids)
1369  _halo_ids[grain._var_index][entity] = grain._var_index;
1370 
1371  for (auto entity : grain._ghosted_ids)
1372  _ghosted_entity_ids[entity] = 1;
1373  }
1374 
1376 }
1377 
1378 void
1380 {
1381  if (_compute_halo_maps)
1382  {
1383  // rank var_index entity_id
1384  std::vector<std::pair<std::size_t, dof_id_type>> halo_ids_all;
1385 
1386  std::vector<int> counts;
1387  std::vector<std::pair<std::size_t, dof_id_type>> local_halo_ids;
1388  std::size_t counter = 0;
1389 
1390  if (_is_master)
1391  {
1392  std::vector<std::vector<std::pair<std::size_t, dof_id_type>>> root_halo_ids(_n_procs);
1393  counts.resize(_n_procs);
1394 
1395  auto & mesh = _mesh.getMesh();
1396  // Loop over the _halo_ids "field" and build minimal lists for all of the other ranks
1397  for (auto var_index = beginIndex(_halo_ids); var_index < _halo_ids.size(); ++var_index)
1398  {
1399  for (const auto & entity_pair : _halo_ids[var_index])
1400  {
1401  DofObject * halo_entity;
1402  if (_is_elemental)
1403  halo_entity = mesh.elem(entity_pair.first);
1404  else
1405  halo_entity = &mesh.node(entity_pair.first);
1406 
1407  root_halo_ids[halo_entity->processor_id()].push_back(
1408  std::make_pair(var_index, entity_pair.first));
1409  }
1410  }
1411 
1412  // Build up the counts vector for MPI scatter
1413  std::size_t global_count = 0;
1414  for (const auto & vector_ref : root_halo_ids)
1415  {
1416  std::copy(vector_ref.begin(), vector_ref.end(), std::back_inserter(halo_ids_all));
1417  counts[counter] = vector_ref.size();
1418  global_count += counts[counter++];
1419  }
1420  }
1421 
1422  _communicator.scatter(halo_ids_all, counts, local_halo_ids);
1423 
1424  // Now add the contributions from the root process to the processor local maps
1425  for (const auto & halo_pair : local_halo_ids)
1426  _halo_ids[halo_pair.first].emplace(std::make_pair(halo_pair.second, halo_pair.first));
1427 
1428  // Finally remove halo markings from stitch regions
1429  for (const auto & grain : _feature_sets)
1430  for (auto local_id : grain._local_ids)
1431  _halo_ids[grain._var_index].erase(local_id);
1432  }
1433 }
1434 
1435 Real
1436 GrainTracker::centroidRegionDistance(std::vector<MeshTools::BoundingBox> & bboxes1,
1437  std::vector<MeshTools::BoundingBox> & bboxes2) const
1438 {
1442  auto min_distance = std::numeric_limits<Real>::max();
1443  for (const auto & bbox1 : bboxes1)
1444  {
1445  const auto centroid_point1 = (bbox1.max() + bbox1.min()) / 2.0;
1446 
1447  for (const auto & bbox2 : bboxes2)
1448  {
1449  const auto centroid_point2 = (bbox2.max() + bbox2.min()) / 2.0;
1450 
1451  // Here we'll calculate a distance between the centroids
1452  auto curr_distance = _mesh.minPeriodicDistance(_var_number, centroid_point1, centroid_point2);
1453 
1454  if (curr_distance < min_distance)
1455  min_distance = curr_distance;
1456  }
1457  }
1458 
1459  return min_distance;
1460 }
1461 
1462 Real
1463 GrainTracker::boundingRegionDistance(std::vector<MeshTools::BoundingBox> & bboxes1,
1464  std::vector<MeshTools::BoundingBox> & bboxes2) const
1465 {
1472  auto min_distance = std::numeric_limits<Real>::max();
1473  for (const auto & bbox1 : bboxes1)
1474  {
1475  for (const auto & bbox2 : bboxes2)
1476  {
1477  // AABB squared distance
1478  Real curr_distance = 0.0;
1479  bool boxes_overlap = true;
1480  for (unsigned int dim = 0; dim < LIBMESH_DIM; ++dim)
1481  {
1482  const auto & min1 = bbox1.min()(dim);
1483  const auto & max1 = bbox1.max()(dim);
1484  const auto & min2 = bbox2.min()(dim);
1485  const auto & max2 = bbox2.max()(dim);
1486 
1487  if (min1 > max2)
1488  {
1489  const auto delta = max2 - min1;
1490  curr_distance += delta * delta;
1491  boxes_overlap = false;
1492  }
1493  else if (min2 > max1)
1494  {
1495  const auto delta = max1 - min2;
1496  curr_distance += delta * delta;
1497  boxes_overlap = false;
1498  }
1499  }
1500 
1501  if (boxes_overlap)
1502  return -1.0; /* all overlaps are treated the same */
1503 
1504  if (curr_distance < min_distance)
1505  min_distance = curr_distance;
1506  }
1507  }
1508 
1509  return min_distance;
1510 }
1511 
1512 unsigned int
1514 {
1522  _max_curr_grain_id = std::max(_max_curr_grain_id + 1,
1523  _reserve_grain_first_index + _n_reserve_ops /* no +1 here!*/);
1524 
1525  return _max_curr_grain_id;
1526 }
1527 
1528 /*************************************************
1529  ************** Helper Structures ****************
1530  ************************************************/
1531 GrainDistance::GrainDistance(Real distance, std::size_t var_index)
1532  : GrainDistance(distance,
1533  var_index,
1534  std::numeric_limits<std::size_t>::max(),
1535  std::numeric_limits<unsigned int>::max())
1536 {
1537 }
1538 
1540  std::size_t var_index,
1541  std::size_t grain_index,
1542  unsigned int grain_id)
1543  : _distance(distance), _var_index(var_index), _grain_index(grain_index), _grain_id(grain_id)
1544 {
1545 }
1546 
1547 bool
1549 {
1550  return _distance < rhs._distance;
1551 }
void remapGrains()
This method is called after trackGrains to remap grains that are too close to each other...
Definition: GrainTracker.C:715
virtual std::size_t getNumberActiveGrains() const override
Returns the number of active grains current stored in the GrainTracker.
Definition: GrainTracker.C:107
virtual Real getEntityValue(dof_id_type entity_id, FieldType field_type, std::size_t var_index=0) const
virtual const std::vector< unsigned int > & getVarToFeatureVector(dof_id_type elem_id) const override
Returns a list of active unique feature ids for a particular element.
Definition: GrainTracker.C:95
void expandEdgeHalos(unsigned int num_layers_to_expand)
This method expands the existing halo set by some width determined by the passed in value...
This class defines the interface for the GrainTracking objects.
const std::size_t _n_vars
static const std::size_t invalid_size_t
Real centroidRegionDistance(std::vector< MeshTools::BoundingBox > &bboxes1, std::vector< MeshTools::BoundingBox > &bboxes2) const
This method returns the minimum periodic distance between the centroids of two vectors of bounding bo...
virtual void execute() override
Definition: GrainTracker.C:173
const std::vector< FeatureData > & getFeatures() const
Return a constant reference to the vector of all discovered features.
const bool _condense_map_info
virtual std::size_t getTotalFeatureCount() const override
Returns the total feature count (active and inactive ids, useful for sizing vectors) ...
Definition: GrainTracker.C:114
static const unsigned int _max_renumbering_recursion
Depth of renumbering recursion (a depth of zero means no recursion)
Definition: GrainTracker.h:183
Status
This enumeration is used to indicate status of the grains in the _unique_grains data structure...
void trackGrains()
On subsequent time_steps, incoming FeatureData objects are compared to previous time_step information...
Definition: GrainTracker.C:398
This object provides the base capability for creating proper polycrystal ICs.
std::map< dof_id_type, std::vector< unsigned int > > _entity_var_to_features
void communicateAndMerge()
This routine handles all of the serialization, communication and deserialization of the data structur...
InputParameters validParams< FeatureFloodCount >()
virtual void initialize() override
unsigned int _old_max_grain_id
The previous max grain id (needed to figure out which ids are new in a given step) ...
Definition: GrainTracker.h:228
void swapSolutionValuesHelper(Node *curr_node, std::size_t curr_var_index, std::size_t new_var_index, std::vector< std::map< Node *, CacheValues >> &cache, RemapCacheMode cache_mode)
Helper method for actually performing the swaps.
virtual void finalize() override
Definition: GrainTracker.C:232
virtual unsigned int getFeatureVar(unsigned int feature_id) const
Returns the variable representing the passed in feature.
std::size_t _var_index
Definition: GrainTracker.h:261
virtual Real getThreshold(std::size_t current_index) const override
Return the starting comparison threshold to use when inspecting an entity during the flood stage...
Definition: GrainTracker.C:188
std::map< dof_id_type, int > _ghosted_entity_ids
The map for holding reconstructed ghosted element information.
This struct is used to hold distance information to other grains in the simulation.
Definition: GrainTracker.h:241
virtual bool doesFeatureIntersectBoundary(unsigned int feature_id) const override
Returns a Boolean indicating whether this feature intersects any boundary.
Definition: GrainTracker.C:139
NonlinearSystemBase & _nl
A reference to the nonlinear system (used for retrieving solution vectors)
Definition: GrainTracker.h:199
void dataLoad(std::istream &stream, GrainTracker::PartialFeatureData &feature, void *context)
Definition: GrainTracker.C:36
void sortAndLabel()
Sort and assign ids to features based on their position in the container after sorting.
virtual void updateFieldInfo() override
This method is used to populate any of the data structures used for storing field data (nodal or elem...
std::size_t _grain_index
Definition: GrainTracker.h:262
void communicateHaloMap()
unsigned int _max_curr_grain_id
Holds the next "regular" grain ID (a grain found or remapped to the standard op vars) ...
Definition: GrainTracker.h:231
GrainTracker(const InputParameters &parameters)
Definition: GrainTracker.C:56
std::vector< std::map< dof_id_type, int > > _halo_ids
The data structure for looking up halos around features.
unsigned int _grain_id
Definition: GrainTracker.h:263
std::vector< MooseVariable * > _vars
The vector of coupled in variables.
bool _first_time
Boolean to indicate the first time this object executes.
Definition: GrainTracker.h:214
virtual void execute() override
bool _is_elemental
Determines if the flood counter is elements or not (nodes)
bool _is_master
Convenience variable for testing master rank.
const PolycrystalUserObjectBase * _poly_ic_uo
An optional IC UserObject which can provide initial data structures to this object.
Definition: GrainTracker.h:208
void prepopulateState(const FeatureFloodCount &ffc_object)
This method extracts the necessary state from the passed in object necessary to continue tracking gra...
Definition: GrainTracker.C:200
void broadcastAndUpdateGrainData()
Broadcast essential Grain information to all processors.
Definition: GrainTracker.C:297
std::vector< FeatureData > _feature_sets
The data structure used to hold the globally unique features.
GrainDistance(Real distance, std::size_t var_index)
unsigned long _var_number
This variable is used to build the periodic node map.
This object will mark nodes or elements of continuous regions all with a unique number for the purpos...
const Real _reserve_op_threshold
The threshold above (or below) where a grain may be found on a reserve op field.
Definition: GrainTracker.h:193
const std::size_t _maps_size
Convenience variable holding the size of all the datastructures size by the number of maps...
virtual void initialize() override
Definition: GrainTracker.C:155
void computeMinDistancesFromGrain(FeatureData &grain, std::vector< std::list< GrainDistance >> &min_distances)
Populates and sorts a min_distances vector with the minimum distances to all grains in the simulation...
Definition: GrainTracker.C:954
void swapSolutionValues(FeatureData &grain, std::size_t new_var_index, std::vector< std::map< Node *, CacheValues >> &cache, RemapCacheMode cache_mode)
A routine for moving all of the solution values from a given grain to a new variable number...
static const unsigned int invalid_id
InputParameters validParams< GrainTrackerInterface >()
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...
InputParameters validParams< GrainTracker >()
Definition: GrainTracker.C:46
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 unsigned int getFeatureVar(unsigned int feature_id) const override
Returns the variable representing the passed in feature.
Definition: GrainTracker.C:101
Real boundingRegionDistance(std::vector< MeshTools::BoundingBox > &bboxes1, std::vector< MeshTools::BoundingBox > &bboxes2) const
This method returns the minimum periodic distance between two vectors of bounding boxes...
virtual Point getGrainCentroid(unsigned int grain_id) const override
Returns the centroid for the given grain number.
Definition: GrainTracker.C:121
std::vector< std::list< FeatureData > > _partial_feature_sets
The data structure used to hold partial and communicated feature data.
const unsigned int _halo_level
The thickness of the halo surrounding each grain.
Definition: GrainTracker.h:180
const bool _compute_halo_maps
Indicates whether or not to communicate halo map information with all ranks.
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.
const bool _remap
Inidicates whether remapping should be done or not (remapping is independent of tracking) ...
Definition: GrainTracker.h:196
virtual Real getVariableValue(unsigned int op_index, const Point &p) const =0
Returns the variable value for a given op_index and mesh point.
bool attemptGrainRenumber(FeatureData &grain, unsigned int depth, unsigned int max_depth)
This is the recursive part of the remapping algorithm.
virtual Real getEntityValue(dof_id_type node_id, FieldType field_type, std::size_t var_index=0) const override
Definition: GrainTracker.C:84
const processor_id_type _n_procs
Convenience variable holding the number of processors in this simulation.
virtual ~GrainTracker()
Definition: GrainTracker.C:81
const bool _is_transient
Boolean to indicate whether this is a Steady or Transient solve.
Definition: GrainTracker.h:234
bool _error_on_grain_creation
Boolean to terminate with an error if a new grain is created during the simulation.
Definition: GrainTracker.h:221
std::vector< unsigned int > _feature_counts_per_map
The number of features seen by this object per map.
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 dataStore(std::ostream &stream, GrainTracker::PartialFeatureData &feature, void *context)
Definition: GrainTracker.C:26
void assignGrains()
When the tracking phase starts (_t_step == _tracking_step) it assigns a unique id to every FeatureDat...
Definition: GrainTracker.C:357
const int _tracking_step
The timestep to begin tracking grains.
Definition: GrainTracker.h:177
unsigned int _reserve_grain_first_index
Holds the first unique grain index when using _reserve_op (all the remaining indices are sequential) ...
Definition: GrainTracker.h:225
const std::size_t _reserve_op_index
The cutoff index where if variable index >= this number, no remapping TO that variable will occur...
Definition: GrainTracker.h:190
MooseMesh & _mesh
A reference to the mesh.
void buildFeatureIdToLocalIndices(unsigned int max_id)
This method builds a lookup map for retrieving the right local feature (by index) given a global inde...
void scatterAndUpdateRanks()
Calls buildLocalToGlobalIndices to build the individual local to global indicies for each rank and sc...
std::vector< FeatureData > & _feature_sets_old
This data structure holds the map of unique grains from the previous time step.
Definition: GrainTracker.h:205
virtual void newGrainCreated(unsigned int new_grain_id)
This method is called when a new grain is detected.
Definition: GrainTracker.C:682
const bool _compute_var_to_feature_map
Indicates whether or not the var to feature map is populated.
static unsigned int counter
const bool _var_index_mode
This variable is used to indicate whether the maps will contain unique region information or just the...
unsigned int getNextUniqueID()
Retrieve the next unique grain number if a new grain is detected during trackGrains.
std::vector< std::map< dof_id_type, int > > _var_index_maps
This map keeps track of which variables own which nodes.
const unsigned short _n_reserve_ops
The number of reserved order parameters.
Definition: GrainTracker.h:186
virtual std::vector< unsigned int > getNewGrainIDs() const override
This method returns all of the new ids generated in an invocation of the GrainTracker.
Definition: GrainTracker.C:703
bool operator<(const GrainDistance &rhs) const