www.mooseframework.org
MultiAppPostprocessorInterpolationTransfer.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 "FEProblem.h"
14 #include "MooseMesh.h"
15 #include "MooseTypes.h"
16 #include "MultiApp.h"
17 
18 #include "libmesh/meshfree_interpolation.h"
19 #include "libmesh/numeric_vector.h"
20 #include "libmesh/system.h"
21 #include "libmesh/radial_basis_interpolation.h"
22 
24 
25 template <>
28 {
30  params.addRequiredParam<AuxVariableName>(
31  "variable", "The auxiliary variable to store the transferred values in.");
32  params.addRequiredParam<PostprocessorName>("postprocessor", "The Postprocessor to interpolate.");
33  params.addParam<unsigned int>(
34  "num_points", 3, "The number of nearest points to use for interpolation.");
35  params.addParam<Real>(
36  "power", 2, "The polynomial power to use for calculation of the decay in the interpolation.");
37 
38  MooseEnum interp_type("inverse_distance radial_basis", "inverse_distance");
39 
40  params.addParam<MooseEnum>("interp_type", interp_type, "The algorithm to use for interpolation.");
41 
42  params.addParam<Real>("radius",
43  -1,
44  "Radius to use for radial_basis interpolation. If negative "
45  "then the radius is taken as the max distance between "
46  "points.");
47 
48  return params;
49 }
50 
52  const InputParameters & parameters)
53  : MultiAppTransfer(parameters),
54  _postprocessor(getParam<PostprocessorName>("postprocessor")),
55  _to_var_name(getParam<AuxVariableName>("variable")),
56  _num_points(getParam<unsigned int>("num_points")),
57  _power(getParam<Real>("power")),
58  _interp_type(getParam<MooseEnum>("interp_type")),
59  _radius(getParam<Real>("radius"))
60 {
61 }
62 
63 void
65 {
66  _console << "Beginning PostprocessorInterpolationTransfer " << name() << std::endl;
67 
68  switch (_direction)
69  {
70  case TO_MULTIAPP:
71  {
72  mooseError("Can't interpolate to a MultiApp!!");
73  break;
74  }
75  case FROM_MULTIAPP:
76  {
77  InverseDistanceInterpolation<LIBMESH_DIM> * idi;
78 
79  switch (_interp_type)
80  {
81  case 0:
82  idi = new InverseDistanceInterpolation<LIBMESH_DIM>(_communicator, _num_points, _power);
83  break;
84  case 1:
85  idi = new RadialBasisInterpolation<LIBMESH_DIM>(_communicator, _radius);
86  break;
87  default:
88  mooseError("Unknown interpolation type!");
89  }
90 
91  std::vector<Point> & src_pts(idi->get_source_points());
92  std::vector<Number> & src_vals(idi->get_source_vals());
93 
94  std::vector<std::string> field_vars;
95  field_vars.push_back(_to_var_name);
96  idi->set_field_variables(field_vars);
97 
98  {
99  for (unsigned int i = 0; i < _multi_app->numGlobalApps(); i++)
100  {
101  if (_multi_app->hasLocalApp(i) && _multi_app->isRootProcessor())
102  {
103  src_pts.push_back(_multi_app->position(i));
104  src_vals.push_back(_multi_app->appPostprocessorValue(i, _postprocessor));
105  }
106  }
107  }
108 
109  // We have only set local values - prepare for use by gathering remote gata
110  idi->prepare_for_use();
111 
112  // Loop over the master nodes and set the value of the variable
113  {
114  System * to_sys = find_sys(_multi_app->problemBase().es(), _to_var_name);
115 
116  unsigned int sys_num = to_sys->number();
117  unsigned int var_num = to_sys->variable_number(_to_var_name);
118 
119  NumericVector<Real> & solution = *to_sys->solution;
120 
121  MooseMesh & mesh = _multi_app->problemBase().mesh();
122 
123  std::vector<std::string> vars;
124 
125  vars.push_back(_to_var_name);
126 
127  for (const auto & node : as_range(mesh.localNodesBegin(), mesh.localNodesEnd()))
128  {
129  if (node->n_dofs(sys_num, var_num) > 0) // If this variable has dofs at this node
130  {
131  std::vector<Point> pts;
132  std::vector<Number> vals;
133 
134  pts.push_back(*node);
135  vals.resize(1);
136 
137  idi->interpolate_field_data(vars, pts, vals);
138 
139  Real value = vals.front();
140 
141  // The zero only works for LAGRANGE!
142  dof_id_type dof = node->dof_number(sys_num, var_num, 0);
143 
144  solution.set(dof, value);
145  }
146  }
147 
148  solution.close();
149  }
150 
151  _multi_app->problemBase().es().update();
152 
153  delete idi;
154 
155  break;
156  }
157  }
158 
159  _console << "Finished PostprocessorInterpolationTransfer " << name() << std::endl;
160 }
InputParameters validParams< MultiAppPostprocessorInterpolationTransfer >()
registerMooseObject("MooseApp", MultiAppPostprocessorInterpolationTransfer)
MeshBase::const_node_iterator localNodesBegin()
Calls local_nodes_begin/end() on the underlying libMesh mesh object.
Definition: MooseMesh.C:2205
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
void mooseError(Args &&... args) const
Definition: MooseObject.h:147
std::shared_ptr< MultiApp > _multi_app
The MultiApp this Transfer is transferring data to or from.
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...
Transfers from spatially varying PostprocessorInterpolations in a MultiApp to the "master" system...
MeshBase::const_node_iterator localNodesEnd()
Definition: MooseMesh.C:2211
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
Definition: MooseMesh.h:74
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:31
InputParameters validParams< MultiAppTransfer >()
Base class for all MultiAppTransfer objects.
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
const MooseEnum _direction
Whether we&#39;re transferring to or from the MultiApp.
MultiAppPostprocessorInterpolationTransfer(const InputParameters &parameters)