www.mooseframework.org
MultiAppMeshFunctionTransfer.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
11 
12 // MOOSE includes
13 #include "DisplacedProblem.h"
14 #include "FEProblem.h"
15 #include "MooseMesh.h"
16 #include "MooseTypes.h"
17 #include "MooseVariableFE.h"
18 
19 #include "libmesh/meshfree_interpolation.h"
20 #include "libmesh/system.h"
21 #include "libmesh/mesh_function.h"
22 #include "libmesh/mesh_tools.h"
23 #include "libmesh/parallel_algebra.h" // for communicator send and receive stuff
24 
26 
27 template <>
30 {
32 
33  params.addParam<bool>(
34  "error_on_miss",
35  false,
36  "Whether or not to error in the case that a target point is not found in the source domain.");
37  return params;
38 }
39 
41  : MultiAppFieldTransfer(parameters), _error_on_miss(getParam<bool>("error_on_miss"))
42 {
43  if (_to_var_names.size() == _from_var_names.size())
44  _var_size = _to_var_names.size();
45  else
46  mooseError("The number of variables to transfer to and from should be equal");
47 }
48 
49 void
51 {
52  _console << "Beginning MeshFunctionTransfer " << name() << std::endl;
53 
54  getAppInfo();
55 
56  _send_points.resize(_var_size);
57  _send_evals.resize(_var_size);
58  _send_ids.resize(_var_size);
59  // loop over the vector of variables and make the transfer one by one
60  for (unsigned int i = 0; i < _var_size; ++i)
62 
63  // Make sure all our sends succeeded.
64  for (unsigned int i = 0; i < _var_size; ++i)
65  for (processor_id_type i_proc = 0; i_proc < n_processors(); ++i_proc)
66  {
67  if (i_proc == processor_id())
68  continue;
69  _send_points[i][i_proc].wait();
70  _send_evals[i][i_proc].wait();
72  _send_ids[i][i_proc].wait();
73  }
74 
75  _console << "Finished MeshFunctionTransfer " << name() << std::endl;
76 
77  postExecute();
78 }
79 
80 void
82 {
83  mooseAssert(i < _var_size, "The variable of index " << i << " does not exist");
84 
93  // Get the bounding boxes for the "from" domains.
94  std::vector<BoundingBox> bboxes = getFromBoundingBoxes();
95 
96  // Figure out how many "from" domains each processor owns.
97  std::vector<unsigned int> froms_per_proc = getFromsPerProc();
98 
99  std::vector<std::vector<Point>> outgoing_points(n_processors());
100  std::vector<std::map<std::pair<unsigned int, unsigned int>, unsigned int>> point_index_map(
101  n_processors());
102  // point_index_map[i_to, element_id] = index
103  // outgoing_points[index] is the first quadrature point in element
104 
105  for (unsigned int i_to = 0; i_to < _to_problems.size(); ++i_to)
106  {
107  System * to_sys = find_sys(*_to_es[i_to], _to_var_names[i]);
108  unsigned int sys_num = to_sys->number();
109  unsigned int var_num = to_sys->variable_number(_to_var_names[i]);
110  MeshBase * to_mesh = &_to_meshes[i_to]->getMesh();
111  bool is_nodal = to_sys->variable_type(var_num).family == LAGRANGE;
112 
113  if (is_nodal)
114  {
115  for (const auto & node : to_mesh->local_node_ptr_range())
116  {
117  // Skip this node if the variable has no dofs at it.
118  if (node->n_dofs(sys_num, var_num) < 1)
119  continue;
120 
121  // Loop over the "froms" on processor i_proc. If the node is found in
122  // any of the "froms", add that node to the vector that will be sent to
123  // i_proc.
124  unsigned int from0 = 0;
125  for (processor_id_type i_proc = 0; i_proc < n_processors();
126  from0 += froms_per_proc[i_proc], ++i_proc)
127  {
128  bool point_found = false;
129  for (unsigned int i_from = from0; i_from < from0 + froms_per_proc[i_proc] && !point_found;
130  ++i_from)
131  {
132  if (bboxes[i_from].contains_point(*node + _to_positions[i_to]))
133  {
134  std::pair<unsigned int, unsigned int> key(i_to, node->id());
135  point_index_map[i_proc][key] = outgoing_points[i_proc].size();
136  outgoing_points[i_proc].push_back(*node + _to_positions[i_to]);
137  point_found = true;
138  }
139  }
140  }
141  }
142  }
143  else // Elemental
144  {
145  for (auto & elem : as_range(to_mesh->local_elements_begin(), to_mesh->local_elements_end()))
146  {
147  Point centroid = elem->centroid();
148 
149  // Skip this element if the variable has no dofs at it.
150  if (elem->n_dofs(sys_num, var_num) < 1)
151  continue;
152 
153  // Loop over the "froms" on processor i_proc. If the elem is found in
154  // any of the "froms", add that elem to the vector that will be sent to
155  // i_proc.
156  unsigned int from0 = 0;
157  for (processor_id_type i_proc = 0; i_proc < n_processors();
158  from0 += froms_per_proc[i_proc], ++i_proc)
159  {
160  bool point_found = false;
161  for (unsigned int i_from = from0; i_from < from0 + froms_per_proc[i_proc] && !point_found;
162  ++i_from)
163  {
164  if (bboxes[i_from].contains_point(centroid + _to_positions[i_to]))
165  {
166  std::pair<unsigned int, unsigned int> key(i_to, elem->id());
167  point_index_map[i_proc][key] = outgoing_points[i_proc].size();
168  outgoing_points[i_proc].push_back(centroid + _to_positions[i_to]);
169  point_found = true;
170  }
171  }
172  }
173  }
174  }
175  }
176 
182  // Get the local bounding boxes.
183  std::vector<BoundingBox> local_bboxes(froms_per_proc[processor_id()]);
184  {
185  // Find the index to the first of this processor's local bounding boxes.
186  unsigned int local_start = 0;
187  for (processor_id_type i_proc = 0; i_proc < n_processors() && i_proc != processor_id();
188  ++i_proc)
189  {
190  local_start += froms_per_proc[i_proc];
191  }
192 
193  // Extract the local bounding boxes.
194  for (unsigned int i_from = 0; i_from < froms_per_proc[processor_id()]; ++i_from)
195  {
196  local_bboxes[i_from] = bboxes[local_start + i_from];
197  }
198  }
199 
200  // Setup the local mesh functions.
201  std::vector<std::shared_ptr<MeshFunction>> local_meshfuns;
202  for (unsigned int i_from = 0; i_from < _from_problems.size(); ++i_from)
203  {
204  FEProblemBase & from_problem = *_from_problems[i_from];
205  MooseVariableFEBase & from_var =
206  from_problem.getVariable(0,
207  _from_var_names[i],
210  System & from_sys = from_var.sys().system();
211  unsigned int from_var_num = from_sys.variable_number(from_var.name());
212 
213  std::shared_ptr<MeshFunction> from_func;
214  // TODO: make MultiAppTransfer give me the right es
215  if (_displaced_source_mesh && from_problem.getDisplacedProblem())
216  from_func.reset(new MeshFunction(from_problem.getDisplacedProblem()->es(),
217  *from_sys.current_local_solution,
218  from_sys.get_dof_map(),
219  from_var_num));
220  else
221  from_func.reset(new MeshFunction(from_problem.es(),
222  *from_sys.current_local_solution,
223  from_sys.get_dof_map(),
224  from_var_num));
225  from_func->init(Trees::ELEMENTS);
226  from_func->enable_out_of_mesh_mode(OutOfMeshValue);
227  local_meshfuns.push_back(from_func);
228  }
229 
230  // Send points to other processors.
231  std::vector<std::vector<Real>> incoming_evals(n_processors());
232  std::vector<std::vector<unsigned int>> incoming_app_ids(n_processors());
233  _send_points[i].resize(n_processors());
234  for (processor_id_type i_proc = 0; i_proc < n_processors(); ++i_proc)
235  {
236  if (i_proc == processor_id())
237  continue;
238  _communicator.send(i_proc, outgoing_points[i_proc], _send_points[i][i_proc]);
239  }
240 
241  // Receive points from other processors, evaluate mesh functions at those
242  // points, and send the values back.
243  _send_evals[i].resize(n_processors());
244  _send_ids[i].resize(n_processors());
245 
246  // Create these here so that they live the entire life of this function
247  // and are NOT reused per processor.
248  std::vector<std::vector<Real>> processor_outgoing_evals(n_processors());
249 
250  for (processor_id_type i_proc = 0; i_proc < n_processors(); ++i_proc)
251  {
252  std::vector<Point> incoming_points;
253  if (i_proc == processor_id())
254  incoming_points = outgoing_points[i_proc];
255  else
256  _communicator.receive(i_proc, incoming_points);
257 
258  std::vector<Real> & outgoing_evals = processor_outgoing_evals[i_proc];
259  outgoing_evals.resize(incoming_points.size(), OutOfMeshValue);
260 
261  std::vector<unsigned int> outgoing_ids(incoming_points.size(), -1); // -1 = largest unsigned int
262  for (unsigned int i_pt = 0; i_pt < incoming_points.size(); ++i_pt)
263  {
264  Point pt = incoming_points[i_pt];
265 
266  // Loop until we've found the lowest-ranked app that actually contains
267  // the quadrature point.
268  for (unsigned int i_from = 0;
269  i_from < _from_problems.size() && outgoing_evals[i_pt] == OutOfMeshValue;
270  ++i_from)
271  {
272  if (local_bboxes[i_from].contains_point(pt))
273  {
274  outgoing_evals[i_pt] = (*local_meshfuns[i_from])(pt - _from_positions[i_from]);
275  if (_direction == FROM_MULTIAPP)
276  outgoing_ids[i_pt] = _local2global_map[i_from];
277  }
278  }
279  }
280 
281  if (i_proc == processor_id())
282  {
283  incoming_evals[i_proc] = outgoing_evals;
284  if (_direction == FROM_MULTIAPP)
285  incoming_app_ids[i_proc] = outgoing_ids;
286  }
287  else
288  {
289  _communicator.send(i_proc, outgoing_evals, _send_evals[i][i_proc]);
290  if (_direction == FROM_MULTIAPP)
291  _communicator.send(i_proc, outgoing_ids, _send_ids[i][i_proc]);
292  }
293  }
294 
302  for (processor_id_type i_proc = 0; i_proc < n_processors(); ++i_proc)
303  {
304  if (i_proc == processor_id())
305  continue;
306 
307  _communicator.receive(i_proc, incoming_evals[i_proc]);
308  if (_direction == FROM_MULTIAPP)
309  _communicator.receive(i_proc, incoming_app_ids[i_proc]);
310  }
311 
312  for (unsigned int i_to = 0; i_to < _to_problems.size(); ++i_to)
313  {
314  System * to_sys = find_sys(*_to_es[i_to], _to_var_names[i]);
315 
316  unsigned int sys_num = to_sys->number();
317  unsigned int var_num = to_sys->variable_number(_to_var_names[i]);
318 
319  NumericVector<Real> * solution = nullptr;
320  switch (_direction)
321  {
322  case TO_MULTIAPP:
323  solution = &getTransferVector(i_to, _to_var_names[i]);
324  break;
325  case FROM_MULTIAPP:
326  solution = to_sys->solution.get();
327  break;
328  default:
329  mooseError("Unknown direction");
330  }
331 
332  MeshBase * to_mesh = &_to_meshes[i_to]->getMesh();
333 
334  bool is_nodal = to_sys->variable_type(var_num).family == LAGRANGE;
335 
336  if (is_nodal)
337  {
338  for (const auto & node : to_mesh->local_node_ptr_range())
339  {
340  // Skip this node if the variable has no dofs at it.
341  if (node->n_dofs(sys_num, var_num) < 1)
342  continue;
343 
344  unsigned int lowest_app_rank = libMesh::invalid_uint;
345  Real best_val = 0.;
346  bool point_found = false;
347  for (unsigned int i_proc = 0; i_proc < incoming_evals.size(); ++i_proc)
348  {
349  // Skip this proc if the node wasn't in it's bounding boxes.
350  std::pair<unsigned int, unsigned int> key(i_to, node->id());
351  if (point_index_map[i_proc].find(key) == point_index_map[i_proc].end())
352  continue;
353  unsigned int i_pt = point_index_map[i_proc][key];
354 
355  // Ignore this proc if it's app has a higher rank than the
356  // previously found lowest app rank.
357  if (_direction == FROM_MULTIAPP)
358  {
359  if (incoming_app_ids[i_proc][i_pt] >= lowest_app_rank)
360  continue;
361  }
362 
363  // Ignore this proc if the point was actually outside its meshes.
364  if (incoming_evals[i_proc][i_pt] == OutOfMeshValue)
365  continue;
366 
367  best_val = incoming_evals[i_proc][i_pt];
368  point_found = true;
369  }
370 
371  if (_error_on_miss && !point_found)
372  mooseError("Point not found! ", *node + _to_positions[i_to]);
373 
374  dof_id_type dof = node->dof_number(sys_num, var_num, 0);
375  solution->set(dof, best_val);
376  }
377  }
378  else // Elemental
379  {
380  for (auto & elem : as_range(to_mesh->local_elements_begin(), to_mesh->local_elements_end()))
381  {
382  // Skip this element if the variable has no dofs at it.
383  if (elem->n_dofs(sys_num, var_num) < 1)
384  continue;
385 
386  unsigned int lowest_app_rank = libMesh::invalid_uint;
387  Real best_val = 0;
388  bool point_found = false;
389  for (unsigned int i_proc = 0; i_proc < incoming_evals.size(); ++i_proc)
390  {
391  // Skip this proc if the elem wasn't in it's bounding boxes.
392  std::pair<unsigned int, unsigned int> key(i_to, elem->id());
393  if (point_index_map[i_proc].find(key) == point_index_map[i_proc].end())
394  continue;
395  unsigned int i_pt = point_index_map[i_proc][key];
396 
397  // Ignore this proc if it's app has a higher rank than the
398  // previously found lowest app rank.
399  if (_direction == FROM_MULTIAPP)
400  {
401  if (incoming_app_ids[i_proc][i_pt] >= lowest_app_rank)
402  continue;
403  }
404 
405  // Ignore this proc if the point was actually outside its meshes.
406  if (incoming_evals[i_proc][i_pt] == OutOfMeshValue)
407  continue;
408 
409  best_val = incoming_evals[i_proc][i_pt];
410  point_found = true;
411  }
412 
413  if (_error_on_miss && !point_found)
414  mooseError("Point not found! ", elem->centroid() + _to_positions[i_to]);
415 
416  dof_id_type dof = elem->dof_number(sys_num, var_num, 0);
417  solution->set(dof, best_val);
418  }
419  }
420  solution->close();
421  to_sys->update();
422  }
423 }
void transferVariable(unsigned int i)
Performs the transfer for the variable of index i.
NumericVector< Real > & getTransferVector(unsigned int i_local, std::string var_name)
If we are transferring to a multiapp, return the appropriate solution vector.
std::vector< std::vector< Parallel::Request > > _send_ids
To send app ids to other processors.
unsigned int _var_size
The number of variables to transfer.
virtual MooseVariableFEBase & getVariable(THREAD_ID tid, const std::string &var_name, Moose::VarKindType expected_var_type=Moose::VarKindType::VAR_ANY, Moose::VarFieldType expected_var_field_type=Moose::VarFieldType::VAR_FIELD_ANY) override
Returns the variable reference for requested variable which must be of the expected_var_type (Nonline...
MultiAppMeshFunctionTransfer(const InputParameters &parameters)
std::vector< std::vector< Parallel::Request > > _send_points
To send points to other processors.
virtual void postExecute()
Add some extra work if necessary after execute().
std::vector< EquationSystems * > _to_es
std::vector< BoundingBox > getFromBoundingBoxes()
Return the bounding boxes of all the "from" domains, including all the domains not local to this proc...
std::vector< FEProblemBase * > _to_problems
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
std::vector< Point > _from_positions
std::vector< AuxVariableName > _to_var_names
void mooseError(Args &&... args) const
Definition: MooseObject.h:147
InputParameters validParams< MultiAppMeshFunctionTransfer >()
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
virtual EquationSystems & es() override
std::vector< std::vector< Parallel::Request > > _send_evals
To send values to other processors.
virtual std::shared_ptr< DisplacedProblem > getDisplacedProblem()
This serves an interface for MultiAppInterpolationTransfer, MultiAppNearestNodeTransfer and so on...
std::vector< unsigned int > getFromsPerProc()
Return the number of "from" domains that each processor owns.
virtual System & system()=0
Get the reference to the libMesh system.
Transfers a vector of variables.
std::vector< unsigned int > _local2global_map
static const Number OutOfMeshValue
Definition: Transfer.h:75
registerMooseObject("MooseApp", MultiAppMeshFunctionTransfer)
std::vector< Point > _to_positions
const std::string & name() const
Get the variable name.
virtual void execute() override
Execute the transfer.
bool _displaced_source_mesh
True if displaced mesh is used for the source mesh, otherwise false.
const std::string & name() const
Get the name of the object.
Definition: MooseObject.h:59
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an option parameter and a documentation string to the InputParameters object...
const ConsoleStream _console
An instance of helper class to write streams to the Console objects.
static System * find_sys(EquationSystems &es, const std::string &var_name)
Small helper function for finding the system containing the variable.
Definition: Transfer.C:65
InputParameters validParams< MultiAppFieldTransfer >()
const MooseEnum _direction
Whether we&#39;re transferring to or from the MultiApp.
void getAppInfo()
This method will fill information into the convenience member variables (_to_problems, _from_meshes, etc.)
SystemBase & sys()
Get the system this variable is part of.
std::vector< VariableName > _from_var_names
std::vector< FEProblemBase * > _from_problems
std::vector< MooseMesh * > _to_meshes