www.mooseframework.org
MultiAppInterpolationTransfer.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 #include "MultiApp.h"
19 
20 #include "libmesh/meshfree_interpolation.h"
21 #include "libmesh/system.h"
22 #include "libmesh/radial_basis_interpolation.h"
23 
25 
26 template <>
29 {
31  params.addClassDescription(
32  "Transfers the value to the target domain from the nearest node in the source domain.");
33  params.addRequiredParam<AuxVariableName>(
34  "variable", "The auxiliary variable to store the transferred values in.");
35  params.addRequiredParam<VariableName>("source_variable", "The variable to transfer from.");
36 
37  params.addParam<unsigned int>(
38  "num_points", 3, "The number of nearest points to use for interpolation.");
39  params.addParam<Real>(
40  "power", 2, "The polynomial power to use for calculation of the decay in the interpolation.");
41 
42  MooseEnum interp_type("inverse_distance radial_basis", "inverse_distance");
43  params.addParam<MooseEnum>("interp_type", interp_type, "The algorithm to use for interpolation.");
44 
45  params.addParam<Real>("radius",
46  -1,
47  "Radius to use for radial_basis interpolation. If negative "
48  "then the radius is taken as the max distance between "
49  "points.");
50 
51  return params;
52 }
53 
55  : MultiAppTransfer(parameters),
56  _to_var_name(getParam<AuxVariableName>("variable")),
57  _from_var_name(getParam<VariableName>("source_variable")),
58  _num_points(getParam<unsigned int>("num_points")),
59  _power(getParam<Real>("power")),
60  _interp_type(getParam<MooseEnum>("interp_type")),
61  _radius(getParam<Real>("radius"))
62 {
63  // This transfer does not work with DistributedMesh
64  _fe_problem.mesh().errorIfDistributedMesh("MultiAppInterpolationTransfer");
65 }
66 
67 void
69 {
70  if (_direction == TO_MULTIAPP)
72  else
74 }
75 
76 void
78 {
79  _console << "Beginning InterpolationTransfer " << name() << std::endl;
80 
81  switch (_direction)
82  {
83  case TO_MULTIAPP:
84  {
85  FEProblemBase & from_problem = _multi_app->problemBase();
86  MooseVariableFEBase & from_var = from_problem.getVariable(
88 
89  MeshBase * from_mesh = NULL;
90 
91  if (_displaced_source_mesh && from_problem.getDisplacedProblem())
92  from_mesh = &from_problem.getDisplacedProblem()->mesh().getMesh();
93  else
94  from_mesh = &from_problem.mesh().getMesh();
95 
96  SystemBase & from_system_base = from_var.sys();
97  System & from_sys = from_system_base.system();
98 
99  unsigned int from_sys_num = from_sys.number();
100  unsigned int from_var_num = from_sys.variable_number(from_var.name());
101 
102  bool from_is_nodal = from_sys.variable_type(from_var_num).family == LAGRANGE;
103 
104  // EquationSystems & from_es = from_sys.get_equation_systems();
105 
106  NumericVector<Number> & from_solution = *from_sys.solution;
107 
108  InverseDistanceInterpolation<LIBMESH_DIM> * idi;
109 
110  switch (_interp_type)
111  {
112  case 0:
113  idi = new InverseDistanceInterpolation<LIBMESH_DIM>(from_sys.comm(), _num_points, _power);
114  break;
115  case 1:
116  idi = new RadialBasisInterpolation<LIBMESH_DIM>(from_sys.comm(), _radius);
117  break;
118  default:
119  mooseError("Unknown interpolation type!");
120  }
121 
122  std::vector<Point> & src_pts(idi->get_source_points());
123  std::vector<Number> & src_vals(idi->get_source_vals());
124 
125  std::vector<std::string> field_vars;
126  field_vars.push_back(_to_var_name);
127  idi->set_field_variables(field_vars);
128 
129  std::vector<std::string> vars;
130  vars.push_back(_to_var_name);
131 
132  if (from_is_nodal)
133  {
134  for (const auto & from_node : from_mesh->local_node_ptr_range())
135  {
136  // Assuming LAGRANGE!
137  if (from_node->n_comp(from_sys_num, from_var_num) == 0)
138  continue;
139 
140  dof_id_type from_dof = from_node->dof_number(from_sys_num, from_var_num, 0);
141 
142  src_pts.push_back(*from_node);
143  src_vals.push_back(from_solution(from_dof));
144  }
145  }
146  else
147  {
148  for (const auto & from_elem :
149  as_range(from_mesh->local_elements_begin(), from_mesh->local_elements_end()))
150  {
151  // Assuming CONSTANT MONOMIAL
152  dof_id_type from_dof = from_elem->dof_number(from_sys_num, from_var_num, 0);
153 
154  src_pts.push_back(from_elem->centroid());
155  src_vals.push_back(from_solution(from_dof));
156  }
157  }
158 
159  // We have only set local values - prepare for use by gathering remote gata
160  idi->prepare_for_use();
161 
162  for (unsigned int i = 0; i < _multi_app->numGlobalApps(); i++)
163  {
164  if (_multi_app->hasLocalApp(i))
165  {
166  Moose::ScopedCommSwapper swapper(_multi_app->comm());
167 
168  // Loop over the master nodes and set the value of the variable
169  System * to_sys = find_sys(_multi_app->appProblemBase(i).es(), _to_var_name);
170 
171  unsigned int sys_num = to_sys->number();
172  unsigned int var_num = to_sys->variable_number(_to_var_name);
173  NumericVector<Real> & solution = _multi_app->appTransferVector(i, _to_var_name);
174 
175  MeshBase * mesh = NULL;
176 
177  if (_displaced_target_mesh && _multi_app->appProblemBase(i).getDisplacedProblem())
178  mesh = &_multi_app->appProblemBase(i).getDisplacedProblem()->mesh().getMesh();
179  else
180  mesh = &_multi_app->appProblemBase(i).mesh().getMesh();
181 
182  bool is_nodal = to_sys->variable_type(var_num).family == LAGRANGE;
183 
184  if (is_nodal)
185  {
186  for (const auto & node : mesh->local_node_ptr_range())
187  {
188  Point actual_position = *node + _multi_app->position(i);
189 
190  if (node->n_dofs(sys_num, var_num) > 0) // If this variable has dofs at this node
191  {
192  std::vector<Point> pts;
193  std::vector<Number> vals;
194 
195  pts.push_back(actual_position);
196  vals.resize(1);
197 
198  idi->interpolate_field_data(vars, pts, vals);
199 
200  Real value = vals.front();
201 
202  // The zero only works for LAGRANGE!
203  if (node->n_comp(from_sys_num, from_var_num) == 0)
204  continue;
205 
206  dof_id_type dof = node->dof_number(sys_num, var_num, 0);
207 
208  solution.set(dof, value);
209  }
210  }
211  }
212  else // Elemental
213  {
214  for (auto & elem : as_range(mesh->local_elements_begin(), mesh->local_elements_end()))
215  {
216  Point centroid = elem->centroid();
217  Point actual_position = centroid + _multi_app->position(i);
218 
219  if (elem->n_dofs(sys_num, var_num) > 0) // If this variable has dofs at this elem
220  {
221  std::vector<Point> pts;
222  std::vector<Number> vals;
223 
224  pts.push_back(actual_position);
225  vals.resize(1);
226 
227  idi->interpolate_field_data(vars, pts, vals);
228 
229  Real value = vals.front();
230 
231  dof_id_type dof = elem->dof_number(sys_num, var_num, 0);
232 
233  solution.set(dof, value);
234  }
235  }
236  }
237 
238  solution.close();
239  to_sys->update();
240  }
241  }
242 
243  delete idi;
244 
245  break;
246  }
247  case FROM_MULTIAPP:
248  {
249  FEProblemBase & to_problem = _multi_app->problemBase();
250  MooseVariableFEBase & to_var = to_problem.getVariable(
252  SystemBase & to_system_base = to_var.sys();
253 
254  System & to_sys = to_system_base.system();
255 
256  NumericVector<Real> & to_solution = *to_sys.solution;
257 
258  unsigned int to_sys_num = to_sys.number();
259 
260  // Only works with a serialized mesh to transfer to!
261  mooseAssert(to_sys.get_mesh().is_serial(),
262  "MultiAppInterpolationTransfer only works with ReplicatedMesh!");
263 
264  unsigned int to_var_num = to_sys.variable_number(to_var.name());
265 
266  // EquationSystems & to_es = to_sys.get_equation_systems();
267 
268  MeshBase * to_mesh = NULL;
269 
270  if (_displaced_target_mesh && to_problem.getDisplacedProblem())
271  to_mesh = &to_problem.getDisplacedProblem()->mesh().getMesh();
272  else
273  to_mesh = &to_problem.mesh().getMesh();
274 
275  bool is_nodal = to_sys.variable_type(to_var_num).family == LAGRANGE;
276 
277  InverseDistanceInterpolation<LIBMESH_DIM> * idi;
278 
279  switch (_interp_type)
280  {
281  case 0:
282  idi = new InverseDistanceInterpolation<LIBMESH_DIM>(to_sys.comm(), _num_points, _power);
283  break;
284  case 1:
285  idi = new RadialBasisInterpolation<LIBMESH_DIM>(to_sys.comm(), _radius);
286  break;
287  default:
288  mooseError("Unknown interpolation type!");
289  }
290 
291  std::vector<Point> & src_pts(idi->get_source_points());
292  std::vector<Number> & src_vals(idi->get_source_vals());
293 
294  std::vector<std::string> field_vars;
295  field_vars.push_back(_to_var_name);
296  idi->set_field_variables(field_vars);
297 
298  std::vector<std::string> vars;
299  vars.push_back(_to_var_name);
300 
301  for (unsigned int i = 0; i < _multi_app->numGlobalApps(); i++)
302  {
303  if (!_multi_app->hasLocalApp(i))
304  continue;
305 
306  Moose::ScopedCommSwapper swapper(_multi_app->comm());
307 
308  FEProblemBase & from_problem = _multi_app->appProblemBase(i);
309  MooseVariableFEBase & from_var =
310  from_problem.getVariable(0,
314  SystemBase & from_system_base = from_var.sys();
315 
316  System & from_sys = from_system_base.system();
317  unsigned int from_sys_num = from_sys.number();
318 
319  unsigned int from_var_num = from_sys.variable_number(from_var.name());
320 
321  bool from_is_nodal = from_sys.variable_type(from_var_num).family == LAGRANGE;
322 
323  // EquationSystems & from_es = from_sys.get_equation_systems();
324 
325  NumericVector<Number> & from_solution = *from_sys.solution;
326 
327  MeshBase * from_mesh = NULL;
328 
329  if (_displaced_source_mesh && from_problem.getDisplacedProblem())
330  from_mesh = &from_problem.getDisplacedProblem()->mesh().getMesh();
331  else
332  from_mesh = &from_problem.mesh().getMesh();
333 
334  Point app_position = _multi_app->position(i);
335 
336  if (from_is_nodal)
337  {
338  for (const auto & from_node : from_mesh->local_node_ptr_range())
339  {
340  // Assuming LAGRANGE!
341  if (from_node->n_comp(from_sys_num, from_var_num) == 0)
342  continue;
343 
344  dof_id_type from_dof = from_node->dof_number(from_sys_num, from_var_num, 0);
345 
346  src_pts.push_back(*from_node + app_position);
347  src_vals.push_back(from_solution(from_dof));
348  }
349  }
350  else
351  {
352  for (auto & from_element :
353  as_range(from_mesh->local_elements_begin(), from_mesh->local_elements_end()))
354  {
355  // Assuming LAGRANGE!
356  if (from_element->n_comp(from_sys_num, from_var_num) == 0)
357  continue;
358 
359  dof_id_type from_dof = from_element->dof_number(from_sys_num, from_var_num, 0);
360 
361  src_pts.push_back(from_element->centroid() + app_position);
362  src_vals.push_back(from_solution(from_dof));
363  }
364  }
365  }
366 
367  // We have only set local values - prepare for use by gathering remote gata
368  idi->prepare_for_use();
369 
370  // Now do the interpolation to the target system
371  if (is_nodal)
372  {
373  for (auto & node : as_range(to_mesh->local_nodes_begin(), to_mesh->local_nodes_end()))
374  {
375  if (node->n_dofs(to_sys_num, to_var_num) > 0) // If this variable has dofs at this node
376  {
377  std::vector<Point> pts;
378  std::vector<Number> vals;
379 
380  pts.push_back(*node);
381  vals.resize(1);
382 
383  idi->interpolate_field_data(vars, pts, vals);
384 
385  Real value = vals.front();
386 
387  // The zero only works for LAGRANGE!
388  dof_id_type dof = node->dof_number(to_sys_num, to_var_num, 0);
389 
390  to_solution.set(dof, value);
391  }
392  }
393  }
394  else // Elemental
395  {
396  for (auto & elem : as_range(to_mesh->local_elements_begin(), to_mesh->local_elements_end()))
397  {
398  Point centroid = elem->centroid();
399 
400  if (elem->n_dofs(to_sys_num, to_var_num) > 0) // If this variable has dofs at this elem
401  {
402  std::vector<Point> pts;
403  std::vector<Number> vals;
404 
405  pts.push_back(centroid);
406  vals.resize(1);
407 
408  idi->interpolate_field_data(vars, pts, vals);
409 
410  Real value = vals.front();
411 
412  dof_id_type dof = elem->dof_number(to_sys_num, to_var_num, 0);
413 
414  to_solution.set(dof, value);
415  }
416  }
417  }
418 
419  to_solution.close();
420  to_sys.update();
421 
422  delete idi;
423 
424  break;
425  }
426  }
427 
428  _console << "Finished InterpolationTransfer " << name() << std::endl;
429 }
430 
431 Node *
433  Real & distance,
434  const MeshBase::const_node_iterator & nodes_begin,
435  const MeshBase::const_node_iterator & nodes_end)
436 {
437  distance = std::numeric_limits<Real>::max();
438  Node * nearest = NULL;
439 
440  for (auto & node : as_range(nodes_begin, nodes_end))
441  {
442  Real current_distance = (p - *node).norm();
443 
444  if (current_distance < distance)
445  {
446  distance = current_distance;
447  nearest = node;
448  }
449  }
450 
451  return nearest;
452 }
InputParameters validParams< MultiAppInterpolationTransfer >()
MultiAppInterpolationTransfer(const InputParameters &parameters)
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...
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:69
void mooseError(Args &&... args) const
Definition: MooseObject.h:147
Base class for a system (of equations)
Definition: SystemBase.h:92
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...
virtual void execute() override
Execute the transfer.
bool _displaced_target_mesh
True if displaced mesh is used for the target mesh, otherwise false.
void errorIfDistributedMesh(std::string name) const
Generate a unified error message if the underlying libMesh mesh is a DistributedMesh.
Definition: MooseMesh.C:2685
void variableIntegrityCheck(const AuxVariableName &var_name) const
Utility to verify that the vEariable in the destination system exists.
virtual std::shared_ptr< DisplacedProblem > getDisplacedProblem()
MeshBase & getMesh()
Accessor for the underlying libMesh Mesh object.
Definition: MooseMesh.C:2567
Copy the value to the target domain from the nearest node in the source domain.
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:31
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.
const std::string & name() const
Get the variable name.
virtual MooseMesh & mesh() override
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 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...
registerMooseObject("MooseApp", MultiAppInterpolationTransfer)
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
const MooseEnum _direction
Whether we&#39;re transferring to or from the MultiApp.
SystemBase & sys()
Get the system this variable is part of.