www.mooseframework.org
MultiAppInterpolationTransfer.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 
16 
17 // MOOSE includes
18 #include "DisplacedProblem.h"
19 #include "FEProblem.h"
20 #include "MooseMesh.h"
21 #include "MooseTypes.h"
22 #include "MooseVariable.h"
23 #include "MultiApp.h"
24 
25 #include "libmesh/meshfree_interpolation.h"
26 #include "libmesh/system.h"
27 #include "libmesh/radial_basis_interpolation.h"
28 
29 template <>
32 {
34  params.addClassDescription(
35  "Transfers the value to the target domain from the nearest node in the source domain.");
36  params.addRequiredParam<AuxVariableName>(
37  "variable", "The auxiliary variable to store the transferred values in.");
38  params.addRequiredParam<VariableName>("source_variable", "The variable to transfer from.");
39  params.addParam<bool>("displaced_source_mesh",
40  false,
41  "Whether or not to use the displaced mesh for the source mesh.");
42  params.addParam<bool>("displaced_target_mesh",
43  false,
44  "Whether or not to use the displaced mesh for the target mesh.");
45 
46  params.addParam<unsigned int>(
47  "num_points", 3, "The number of nearest points to use for interpolation.");
48  params.addParam<Real>(
49  "power", 2, "The polynomial power to use for calculation of the decay in the interpolation.");
50 
51  MooseEnum interp_type("inverse_distance radial_basis", "inverse_distance");
52  params.addParam<MooseEnum>("interp_type", interp_type, "The algorithm to use for interpolation.");
53 
54  params.addParam<Real>("radius",
55  -1,
56  "Radius to use for radial_basis interpolation. If negative "
57  "then the radius is taken as the max distance between "
58  "points.");
59 
60  return params;
61 }
62 
64  : MultiAppTransfer(parameters),
65  _to_var_name(getParam<AuxVariableName>("variable")),
66  _from_var_name(getParam<VariableName>("source_variable")),
67  _num_points(getParam<unsigned int>("num_points")),
68  _power(getParam<Real>("power")),
69  _interp_type(getParam<MooseEnum>("interp_type")),
70  _radius(getParam<Real>("radius"))
71 {
72  // This transfer does not work with DistributedMesh
73  _fe_problem.mesh().errorIfDistributedMesh("MultiAppInterpolationTransfer");
74  _displaced_source_mesh = getParam<bool>("displaced_source_mesh");
75  _displaced_target_mesh = getParam<bool>("displaced_target_mesh");
76 }
77 
78 void
80 {
81  if (_direction == TO_MULTIAPP)
83  else
85 }
86 
87 void
89 {
90  _console << "Beginning InterpolationTransfer " << name() << std::endl;
91 
92  switch (_direction)
93  {
94  case TO_MULTIAPP:
95  {
96  FEProblemBase & from_problem = _multi_app->problemBase();
97  MooseVariable & from_var = from_problem.getVariable(0, _from_var_name);
98 
99  MeshBase * from_mesh = NULL;
100 
101  if (_displaced_source_mesh && from_problem.getDisplacedProblem())
102  from_mesh = &from_problem.getDisplacedProblem()->mesh().getMesh();
103  else
104  from_mesh = &from_problem.mesh().getMesh();
105 
106  SystemBase & from_system_base = from_var.sys();
107  System & from_sys = from_system_base.system();
108 
109  unsigned int from_sys_num = from_sys.number();
110  unsigned int from_var_num = from_sys.variable_number(from_var.name());
111 
112  bool from_is_nodal = from_sys.variable_type(from_var_num).family == LAGRANGE;
113 
114  // EquationSystems & from_es = from_sys.get_equation_systems();
115 
116  NumericVector<Number> & from_solution = *from_sys.solution;
117 
118  InverseDistanceInterpolation<LIBMESH_DIM> * idi;
119 
120  switch (_interp_type)
121  {
122  case 0:
123  idi = new InverseDistanceInterpolation<LIBMESH_DIM>(from_sys.comm(), _num_points, _power);
124  break;
125  case 1:
126  idi = new RadialBasisInterpolation<LIBMESH_DIM>(from_sys.comm(), _radius);
127  break;
128  default:
129  mooseError("Unknown interpolation type!");
130  }
131 
132  std::vector<Point> & src_pts(idi->get_source_points());
133  std::vector<Number> & src_vals(idi->get_source_vals());
134 
135  std::vector<std::string> field_vars;
136  field_vars.push_back(_to_var_name);
137  idi->set_field_variables(field_vars);
138 
139  std::vector<std::string> vars;
140  vars.push_back(_to_var_name);
141 
142  if (from_is_nodal)
143  {
144  MeshBase::const_node_iterator from_nodes_it = from_mesh->local_nodes_begin();
145  MeshBase::const_node_iterator from_nodes_end = from_mesh->local_nodes_end();
146 
147  for (; from_nodes_it != from_nodes_end; ++from_nodes_it)
148  {
149  Node * from_node = *from_nodes_it;
150 
151  // Assuming LAGRANGE!
152  dof_id_type from_dof = from_node->dof_number(from_sys_num, from_var_num, 0);
153 
154  src_pts.push_back(*from_node);
155  src_vals.push_back(from_solution(from_dof));
156  }
157  }
158  else
159  {
160  MeshBase::const_element_iterator from_elements_it = from_mesh->local_elements_begin();
161  MeshBase::const_element_iterator from_elements_end = from_mesh->local_elements_end();
162 
163  for (; from_elements_it != from_elements_end; ++from_elements_it)
164  {
165  Elem * from_elem = *from_elements_it;
166 
167  // Assuming CONSTANT MONOMIAL
168  dof_id_type from_dof = from_elem->dof_number(from_sys_num, from_var_num, 0);
169 
170  src_pts.push_back(from_elem->centroid());
171  src_vals.push_back(from_solution(from_dof));
172  }
173  }
174 
175  // We have only set local values - prepare for use by gathering remote gata
176  idi->prepare_for_use();
177 
178  for (unsigned int i = 0; i < _multi_app->numGlobalApps(); i++)
179  {
180  if (_multi_app->hasLocalApp(i))
181  {
182  Moose::ScopedCommSwapper swapper(_multi_app->comm());
183 
184  // Loop over the master nodes and set the value of the variable
185  System * to_sys = find_sys(_multi_app->appProblemBase(i).es(), _to_var_name);
186 
187  unsigned int sys_num = to_sys->number();
188  unsigned int var_num = to_sys->variable_number(_to_var_name);
189  NumericVector<Real> & solution = _multi_app->appTransferVector(i, _to_var_name);
190 
191  MeshBase * mesh = NULL;
192 
193  if (_displaced_target_mesh && _multi_app->appProblemBase(i).getDisplacedProblem())
194  mesh = &_multi_app->appProblemBase(i).getDisplacedProblem()->mesh().getMesh();
195  else
196  mesh = &_multi_app->appProblemBase(i).mesh().getMesh();
197 
198  bool is_nodal = to_sys->variable_type(var_num).family == LAGRANGE;
199 
200  if (is_nodal)
201  {
202  MeshBase::const_node_iterator node_it = mesh->local_nodes_begin();
203  MeshBase::const_node_iterator node_end = mesh->local_nodes_end();
204 
205  for (; node_it != node_end; ++node_it)
206  {
207  Node * node = *node_it;
208 
209  Point actual_position = *node + _multi_app->position(i);
210 
211  if (node->n_dofs(sys_num, var_num) > 0) // If this variable has dofs at this node
212  {
213  std::vector<Point> pts;
214  std::vector<Number> vals;
215 
216  pts.push_back(actual_position);
217  vals.resize(1);
218 
219  idi->interpolate_field_data(vars, pts, vals);
220 
221  Real value = vals.front();
222 
223  // The zero only works for LAGRANGE!
224  dof_id_type dof = node->dof_number(sys_num, var_num, 0);
225 
226  solution.set(dof, value);
227  }
228  }
229  }
230  else // Elemental
231  {
232  MeshBase::const_element_iterator elem_it = mesh->local_elements_begin();
233  MeshBase::const_element_iterator elem_end = mesh->local_elements_end();
234 
235  for (; elem_it != elem_end; ++elem_it)
236  {
237  Elem * elem = *elem_it;
238 
239  Point centroid = elem->centroid();
240  Point actual_position = centroid + _multi_app->position(i);
241 
242  if (elem->n_dofs(sys_num, var_num) > 0) // If this variable has dofs at this elem
243  {
244  std::vector<Point> pts;
245  std::vector<Number> vals;
246 
247  pts.push_back(actual_position);
248  vals.resize(1);
249 
250  idi->interpolate_field_data(vars, pts, vals);
251 
252  Real value = vals.front();
253 
254  dof_id_type dof = elem->dof_number(sys_num, var_num, 0);
255 
256  solution.set(dof, value);
257  }
258  }
259  }
260 
261  solution.close();
262  to_sys->update();
263  }
264  }
265 
266  delete idi;
267 
268  break;
269  }
270  case FROM_MULTIAPP:
271  {
272  FEProblemBase & to_problem = _multi_app->problemBase();
273  MooseVariable & to_var = to_problem.getVariable(0, _to_var_name);
274  SystemBase & to_system_base = to_var.sys();
275 
276  System & to_sys = to_system_base.system();
277 
278  NumericVector<Real> & to_solution = *to_sys.solution;
279 
280  unsigned int to_sys_num = to_sys.number();
281 
282  // Only works with a serialized mesh to transfer to!
283  mooseAssert(to_sys.get_mesh().is_serial(),
284  "MultiAppInterpolationTransfer only works with ReplicatedMesh!");
285 
286  unsigned int to_var_num = to_sys.variable_number(to_var.name());
287 
288  // EquationSystems & to_es = to_sys.get_equation_systems();
289 
290  MeshBase * to_mesh = NULL;
291 
292  if (_displaced_target_mesh && to_problem.getDisplacedProblem())
293  to_mesh = &to_problem.getDisplacedProblem()->mesh().getMesh();
294  else
295  to_mesh = &to_problem.mesh().getMesh();
296 
297  bool is_nodal = to_sys.variable_type(to_var_num).family == LAGRANGE;
298 
299  InverseDistanceInterpolation<LIBMESH_DIM> * idi;
300 
301  switch (_interp_type)
302  {
303  case 0:
304  idi = new InverseDistanceInterpolation<LIBMESH_DIM>(to_sys.comm(), _num_points, _power);
305  break;
306  case 1:
307  idi = new RadialBasisInterpolation<LIBMESH_DIM>(to_sys.comm(), _radius);
308  break;
309  default:
310  mooseError("Unknown interpolation type!");
311  }
312 
313  std::vector<Point> & src_pts(idi->get_source_points());
314  std::vector<Number> & src_vals(idi->get_source_vals());
315 
316  std::vector<std::string> field_vars;
317  field_vars.push_back(_to_var_name);
318  idi->set_field_variables(field_vars);
319 
320  std::vector<std::string> vars;
321  vars.push_back(_to_var_name);
322 
323  for (unsigned int i = 0; i < _multi_app->numGlobalApps(); i++)
324  {
325  if (!_multi_app->hasLocalApp(i))
326  continue;
327 
328  Moose::ScopedCommSwapper swapper(_multi_app->comm());
329 
330  FEProblemBase & from_problem = _multi_app->appProblemBase(i);
331  MooseVariable & from_var = from_problem.getVariable(0, _from_var_name);
332  SystemBase & from_system_base = from_var.sys();
333 
334  System & from_sys = from_system_base.system();
335  unsigned int from_sys_num = from_sys.number();
336 
337  unsigned int from_var_num = from_sys.variable_number(from_var.name());
338 
339  bool from_is_nodal = from_sys.variable_type(from_var_num).family == LAGRANGE;
340 
341  // EquationSystems & from_es = from_sys.get_equation_systems();
342 
343  NumericVector<Number> & from_solution = *from_sys.solution;
344 
345  MeshBase * from_mesh = NULL;
346 
347  if (_displaced_source_mesh && from_problem.getDisplacedProblem())
348  from_mesh = &from_problem.getDisplacedProblem()->mesh().getMesh();
349  else
350  from_mesh = &from_problem.mesh().getMesh();
351 
352  Point app_position = _multi_app->position(i);
353 
354  if (from_is_nodal)
355  {
356  MeshBase::const_node_iterator from_nodes_it = from_mesh->local_nodes_begin();
357  MeshBase::const_node_iterator from_nodes_end = from_mesh->local_nodes_end();
358 
359  for (; from_nodes_it != from_nodes_end; ++from_nodes_it)
360  {
361  Node * from_node = *from_nodes_it;
362 
363  // Assuming LAGRANGE!
364  dof_id_type from_dof = from_node->dof_number(from_sys_num, from_var_num, 0);
365 
366  src_pts.push_back(*from_node + app_position);
367  src_vals.push_back(from_solution(from_dof));
368  }
369  }
370  else
371  {
372  MeshBase::const_element_iterator from_elements_it = from_mesh->local_elements_begin();
373  MeshBase::const_element_iterator from_elements_end = from_mesh->local_elements_end();
374 
375  for (; from_elements_it != from_elements_end; ++from_elements_it)
376  {
377  Elem * from_element = *from_elements_it;
378 
379  // Assuming LAGRANGE!
380  dof_id_type from_dof = from_element->dof_number(from_sys_num, from_var_num, 0);
381 
382  src_pts.push_back(from_element->centroid() + app_position);
383  src_vals.push_back(from_solution(from_dof));
384  }
385  }
386  }
387 
388  // We have only set local values - prepare for use by gathering remote gata
389  idi->prepare_for_use();
390 
391  // Now do the interpolation to the target system
392  if (is_nodal)
393  {
394  MeshBase::const_node_iterator node_it = to_mesh->local_nodes_begin();
395  MeshBase::const_node_iterator node_end = to_mesh->local_nodes_end();
396 
397  for (; node_it != node_end; ++node_it)
398  {
399  Node * node = *node_it;
400 
401  if (node->n_dofs(to_sys_num, to_var_num) > 0) // If this variable has dofs at this node
402  {
403  std::vector<Point> pts;
404  std::vector<Number> vals;
405 
406  pts.push_back(*node);
407  vals.resize(1);
408 
409  idi->interpolate_field_data(vars, pts, vals);
410 
411  Real value = vals.front();
412 
413  // The zero only works for LAGRANGE!
414  dof_id_type dof = node->dof_number(to_sys_num, to_var_num, 0);
415 
416  to_solution.set(dof, value);
417  }
418  }
419  }
420  else // Elemental
421  {
422  MeshBase::const_element_iterator elem_it = to_mesh->local_elements_begin();
423  MeshBase::const_element_iterator elem_end = to_mesh->local_elements_end();
424 
425  for (; elem_it != elem_end; ++elem_it)
426  {
427  Elem * elem = *elem_it;
428 
429  Point centroid = elem->centroid();
430 
431  if (elem->n_dofs(to_sys_num, to_var_num) > 0) // If this variable has dofs at this elem
432  {
433  std::vector<Point> pts;
434  std::vector<Number> vals;
435 
436  pts.push_back(centroid);
437  vals.resize(1);
438 
439  idi->interpolate_field_data(vars, pts, vals);
440 
441  Real value = vals.front();
442 
443  dof_id_type dof = elem->dof_number(to_sys_num, to_var_num, 0);
444 
445  to_solution.set(dof, value);
446  }
447  }
448  }
449 
450  to_solution.close();
451  to_sys.update();
452 
453  delete idi;
454 
455  break;
456  }
457  }
458 
459  _console << "Finished InterpolationTransfer " << name() << std::endl;
460 }
461 
462 Node *
464  Real & distance,
465  const MeshBase::const_node_iterator & nodes_begin,
466  const MeshBase::const_node_iterator & nodes_end)
467 {
468  distance = std::numeric_limits<Real>::max();
469  Node * nearest = NULL;
470 
471  for (MeshBase::const_node_iterator node_it = nodes_begin; node_it != nodes_end; ++node_it)
472  {
473  Real current_distance = (p - *(*node_it)).norm();
474 
475  if (current_distance < distance)
476  {
477  distance = current_distance;
478  nearest = *node_it;
479  }
480  }
481 
482  return nearest;
483 }
const std::string & name() const
Get the name of the object.
Definition: MooseObject.h:47
InputParameters validParams< MultiAppInterpolationTransfer >()
void variableIntegrityCheck(const AuxVariableName &var_name) const
Utility to verify that the vEariable in the destination system exists.
MultiAppInterpolationTransfer(const InputParameters &parameters)
Class for stuff related to variables.
Definition: MooseVariable.h:43
Node * getNearestNode(const Point &p, Real &distance, const MeshBase::const_node_iterator &nodes_begin, const MeshBase::const_node_iterator &nodes_end)
Return the nearest node to the point p.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
FEProblemBase & _fe_problem
Definition: Transfer.h:75
void errorIfDistributedMesh(std::string name) const
Generate a unified error message if the underlying libMesh mesh is a DistributedMesh.
Definition: MooseMesh.C:2526
Base class for a system (of equations)
Definition: SystemBase.h:91
std::shared_ptr< MultiApp > _multi_app
The MultiApp this Transfer is transferring data to or from.
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
void addRequiredParam(const std::string &name, const std::string &doc_string)
This method adds a parameter and documentation string to the InputParameters object that will be extr...
const std::string & name() const
Get the variable name.
virtual void execute() override
Execute the transfer.
virtual std::shared_ptr< DisplacedProblem > getDisplacedProblem()
MeshBase & getMesh()
Accessor for the underlying libMesh Mesh object.
Definition: MooseMesh.C:2408
MooseEnum _direction
Whether we&#39;re transferring to or from the MultiApp.
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:37
virtual System & system()=0
Get the reference to the libMesh system.
virtual void initialSetup() override
Method called at the beginning of the simulation for checking integrity or doing one-time setup...
InputParameters validParams< MultiAppTransfer >()
Base class for all MultiAppTransfer objects.
virtual MooseMesh & mesh() override
void addClassDescription(const std::string &doc_string)
This method adds a description of the class that will be displayed in the input file syntax dump...
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...
void mooseError(Args &&...args) const
Definition: MooseObject.h:80
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:70
SystemBase & sys()
Get the system this variable is part of.
virtual MooseVariable & getVariable(THREAD_ID tid, const std::string &var_name) override
Returns the variable reference for requested variable which may be in any system. ...