www.mooseframework.org
DiscreteNucleationInserter.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 
9 #include "libmesh/parallel_algebra.h"
10 
11 #include "libmesh/quadrature.h"
12 
13 template <>
14 InputParameters
16 {
17  InputParameters params = validParams<ElementUserObject>();
18  params.addClassDescription("Manages the list of currently active nucleation sites and adds new "
19  "sites according to a given probability function.");
20  params.addRequiredParam<MaterialPropertyName>(
21  "probability", "Probability density for inserting a discrete nucleus");
22  params.addRequiredParam<Real>("hold_time", "Time to keep each nucleus active");
23  params.addParam<Point>("test", "Insert a fixed nucleus at a point in the simulation cell");
24  MultiMooseEnum setup_options(SetupInterface::getExecuteOptions());
25  setup_options = "timestep_end";
26  params.set<MultiMooseEnum>("execute_on") = setup_options;
27  return params;
28 }
29 
30 DiscreteNucleationInserter::DiscreteNucleationInserter(const InputParameters & parameters)
31  : ElementUserObject(parameters),
32  _probability(getMaterialProperty<Real>("probability")),
33  _hold_time(getParam<Real>("hold_time")),
34  _changes_made(0),
35  _global_nucleus_list(declareRestartableData("global_nucleus_list", NucleusList(0))),
36  _local_nucleus_list(declareRestartableData("local_nucleus_list", NucleusList(0)))
37 {
38  setRandomResetFrequency(EXEC_TIMESTEP_END);
39 
40  // debugging code (this will insert the entry into every processors list, but duplicate entries in
41  // global should be OK)
42  // we also assume that time starts at 0! But hey, this is only for debugging anyways...
43  if (isParamValid("test"))
44  _insert_test = true;
45  else
46  _insert_test = false;
47 
48  // force a map rebuild after restart or recover
49  _changes_made = _app.isRecovering() || _app.isRestarting();
50 }
51 
52 void
54 {
55  _changes_made = 0;
56 
57  // insert test nucleus once
58  if (_insert_test)
59  {
60  _local_nucleus_list.push_back(NucleusLocation(_hold_time, getParam<Point>("test")));
61  _changes_made++;
62  _insert_test = false;
63  }
64 
65  // expire entries from the local nucleus list (if the current timestep converged)
66  if (_fe_problem.converged())
67  {
68  unsigned int i = 0;
69  while (i < _local_nucleus_list.size())
70  {
71  if (_local_nucleus_list[i].first + _fe_problem.dt() <= _fe_problem.time())
72  {
73  // remove entry (by replacing with last element and shrinking size by one)
75  _local_nucleus_list.pop_back();
76  _changes_made++;
77  }
78  else
79  ++i;
80  }
81  }
82 
83  // we reassemble this list at every timestep
84  _global_nucleus_list.clear();
85 }
86 
87 void
89 {
90  // check each qp for potential nucleation
91  // TODO: we might as well place the nuclei at random positions within the element...
92  for (unsigned int qp = 0; qp < _qrule->n_points(); ++qp)
93  if (getRandomReal() < _probability[qp] * _JxW[qp] * _coord[qp] * _fe_problem.dt())
94  {
95  _local_nucleus_list.push_back(
96  NucleusLocation(_fe_problem.dt() + _fe_problem.time() + _hold_time, _q_point[qp]));
97  _changes_made++;
98  }
99 }
100 
101 void
103 {
104  // combine _local_nucleus_list entries from all threads on the current process
105  const DiscreteNucleationInserter & uo = static_cast<const DiscreteNucleationInserter &>(y);
106  _global_nucleus_list.insert(
107  _global_nucleus_list.end(), uo._local_nucleus_list.begin(), uo._local_nucleus_list.end());
109 }
110 
111 void
113 {
114  // add the _local_nucleus_list of thread zero
115  _global_nucleus_list.insert(
117 
123  std::vector<Real> comm_buffer(_global_nucleus_list.size() * 4);
124  for (unsigned i = 0; i < _global_nucleus_list.size(); ++i)
125  {
126  comm_buffer[i * 4 + 0] = _global_nucleus_list[i].first;
127  comm_buffer[i * 4 + 1] = _global_nucleus_list[i].second(0);
128  comm_buffer[i * 4 + 2] = _global_nucleus_list[i].second(1);
129  comm_buffer[i * 4 + 3] = _global_nucleus_list[i].second(2);
130  }
131 
132  // combine _global_nucleus_lists from all MPI ranks
133  _communicator.allgather(comm_buffer);
134 
135  // unpack the gathered _global_nucleus_list
136  unsigned int n = comm_buffer.size() / 4;
137  mooseAssert(comm_buffer.size() % 4 == 0,
138  "Communication buffer has an unexpected size (not divisible by 4)");
139  _global_nucleus_list.resize(n);
140  for (unsigned i = 0; i < n; ++i)
141  {
142  _global_nucleus_list[i].first = comm_buffer[i * 4 + 0];
143  _global_nucleus_list[i].second(0) = comm_buffer[i * 4 + 1];
144  _global_nucleus_list[i].second(1) = comm_buffer[i * 4 + 2];
145  _global_nucleus_list[i].second(2) = comm_buffer[i * 4 + 3];
146  }
147 
148  // get the global number of changes (i.e. changes to _global_nucleus_list)
149  gatherSum(_changes_made);
150 }
const MaterialProperty< Real > & _probability
Nucleation rate density (should be a material property implementing nucleation theory) ...
unsigned int _changes_made
count the number of nucleus deletions and insertions
This UserObject manages the insertion and expiration of nuclei in the simulation domain it manages a ...
DiscreteNucleationInserter(const InputParameters &parameters)
virtual void threadJoin(const UserObject &y)
std::vector< NucleusLocation > NucleusList
Every MPI task should keep a full list of nuclei (in case they cross domains with their finite radii)...
InputParameters validParams< DiscreteNucleationInserter >()
std::pair< Real, Point > NucleusLocation
A nucleus has an expiration time and a location.
NucleusList & _local_nucleus_list
the local nucleus list of nuclei centered in the domain of the current processor
Real _hold_time
Duration of time each nucleus is kept active after insertion.
bool _insert_test
insert test location
NucleusList & _global_nucleus_list
the global list of all nuclei over all processors