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 #include "MooseAppCoordTransform.h"
18 
19 #include "libmesh/meshfree_interpolation.h"
20 #include "libmesh/numeric_vector.h"
21 #include "libmesh/system.h"
22 #include "libmesh/radial_basis_interpolation.h"
23 
25 
28 {
30  params.addClassDescription("Transfer postprocessor data from sub-application into field data on "
31  "the parent application.");
32  params.addRequiredParam<AuxVariableName>(
33  "variable", "The auxiliary variable to store the transferred values in.");
34  params.addRequiredParam<PostprocessorName>("postprocessor", "The Postprocessor to interpolate.");
35  params.addParam<unsigned int>(
36  "num_points", 3, "The number of nearest points to use for interpolation.");
37  params.addParam<Real>(
38  "power", 2, "The polynomial power to use for calculation of the decay in the interpolation.");
39 
40  MooseEnum interp_type("inverse_distance radial_basis", "inverse_distance");
41 
42  params.addParam<MooseEnum>("interp_type", interp_type, "The algorithm to use for interpolation.");
43 
44  params.addParam<Real>("radius",
45  -1,
46  "Radius to use for radial_basis interpolation. If negative "
47  "then the radius is taken as the max distance between "
48  "points.");
49  return params;
50 }
51 
53  const InputParameters & parameters)
54  : MultiAppTransfer(parameters),
55  _postprocessor(getParam<PostprocessorName>("postprocessor")),
56  _to_var_name(getParam<AuxVariableName>("variable")),
57  _num_points(getParam<unsigned int>("num_points")),
58  _power(getParam<Real>("power")),
59  _interp_type(getParam<MooseEnum>("interp_type")),
60  _radius(getParam<Real>("radius")),
61  _nodal(false)
62 {
63  if (isParamValid("to_multi_app"))
64  paramError("to_multi_app", "Unused parameter; only from-MultiApp transfers are implemented");
65 
66  auto & to_fe_type =
67  getFromMultiApp()->problemBase().getStandardVariable(0, _to_var_name).feType();
68  if ((to_fe_type.order != CONSTANT || to_fe_type.family != MONOMIAL) &&
69  (to_fe_type.order != FIRST || to_fe_type.family != LAGRANGE))
70  paramError("variable", "Must be either CONSTANT MONOMIAL or FIRST LAGRANGE");
71 
72  _nodal = to_fe_type.family == LAGRANGE;
73 }
74 
75 void
77 {
78  TIME_SECTION("MultiAppPostprocessorInterpolationTransfer::execute()",
79  5,
80  "Transferring/interpolating postprocessors");
81 
82  switch (_current_direction)
83  {
84  case TO_MULTIAPP:
85  {
86  mooseError("Interpolation from a variable to a MultiApp postprocessors has not been "
87  "implemented. Use MultiAppVariableValueSamplePostprocessorTransfer!");
88  break;
89  }
90  case FROM_MULTIAPP:
91  {
92  std::unique_ptr<InverseDistanceInterpolation<LIBMESH_DIM>> idi;
93 
94  switch (_interp_type)
95  {
96  case 0:
97  idi = std::make_unique<InverseDistanceInterpolation<LIBMESH_DIM>>(
99  break;
100  case 1:
101  idi = std::make_unique<RadialBasisInterpolation<LIBMESH_DIM>>(_communicator, _radius);
102  break;
103  default:
104  mooseError("Unknown interpolation type!");
105  }
106 
107  std::vector<Point> & src_pts(idi->get_source_points());
108  std::vector<Number> & src_vals(idi->get_source_vals());
109 
110  std::vector<std::string> field_vars;
111  field_vars.push_back(_to_var_name);
112  idi->set_field_variables(field_vars);
113 
114  {
115  mooseAssert(_to_transforms.size() == 1, "There should only be one transform here");
116  const auto & to_coord_transform = *_to_transforms[0];
117  for (unsigned int i = 0; i < getFromMultiApp()->numGlobalApps(); i++)
118  {
119  if (getFromMultiApp()->hasLocalApp(i) && getFromMultiApp()->isRootProcessor())
120  {
121  // Evaluation of the _from_transform at the origin yields the transformed position of
122  // the from multi-app
123  if (!getFromMultiApp()->runningInPosition())
124  src_pts.push_back(to_coord_transform.mapBack((*_from_transforms[i])(Point(0))));
125  else
126  // if running in position, the subapp mesh has been transformed so the translation
127  // is no longer applied by the transform
128  src_pts.push_back(to_coord_transform.mapBack(
129  (*_from_transforms[i])(getFromMultiApp()->position(i))));
130  src_vals.push_back(getFromMultiApp()->appPostprocessorValue(i, _postprocessor));
131  }
132  }
133  }
134 
135  // We have only set local values - prepare for use by gathering remote gata
136  idi->prepare_for_use();
137 
138  // Loop over the parent app nodes and set the value of the variable
139  {
140  System * to_sys = find_sys(getFromMultiApp()->problemBase().es(), _to_var_name);
141 
142  unsigned int sys_num = to_sys->number();
143  unsigned int var_num = to_sys->variable_number(_to_var_name);
144 
145  NumericVector<Real> & solution = *to_sys->solution;
146 
147  MooseMesh & mesh = getFromMultiApp()->problemBase().mesh();
148 
149  std::vector<std::string> vars;
150 
151  vars.push_back(_to_var_name);
152 
153  if (_nodal)
154  {
155  // handle linear lagrange shape functions
156  for (const auto & node : as_range(mesh.localNodesBegin(), mesh.localNodesEnd()))
157  {
158  if (node->n_dofs(sys_num, var_num) > 0) // If this variable has dofs at this node
159  {
160  std::vector<Point> pts;
161  std::vector<Number> vals;
162 
163  pts.push_back(*node);
164  vals.resize(1);
165 
166  idi->interpolate_field_data(vars, pts, vals);
167 
168  Real value = vals.front();
169 
170  // The zero only works for LAGRANGE!
171  dof_id_type dof = node->dof_number(sys_num, var_num, 0);
172 
173  solution.set(dof, value);
174  }
175  }
176  }
177  else
178  {
179  // handle constant monomial shape functions
180  for (const auto & elem :
181  as_range(mesh.getMesh().local_elements_begin(), mesh.getMesh().local_elements_end()))
182  {
183  // Exclude the elements without dofs
184  if (elem->n_dofs(sys_num, var_num) > 0)
185  {
186  std::vector<Point> pts;
187  std::vector<Number> vals;
188 
189  pts.push_back(elem->vertex_average());
190  vals.resize(1);
191 
192  idi->interpolate_field_data(vars, pts, vals);
193 
194  Real value = vals.front();
195 
196  dof_id_type dof = elem->dof_number(sys_num, var_num, 0);
197  solution.set(dof, value);
198  }
199  }
200  }
201 
202  solution.close();
203  }
204 
205  getFromMultiApp()->problemBase().es().update();
206 
207  break;
208  }
209  }
210 }
LAGRANGE
std::vector< std::unique_ptr< MultiAppCoordTransform > > _to_transforms
const std::shared_ptr< MultiApp > getFromMultiApp() const
Get the MultiApp to transfer data from.
MooseEnum _current_direction
Definition: Transfer.h:106
registerMooseObject("MooseApp", MultiAppPostprocessorInterpolationTransfer)
FIRST
MeshBase & mesh
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
Real _radius
Radius to consider for inverse interpolation.
const Parallel::Communicator & _communicator
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 postprocessors in child apps of a MultiApp in different locations in parent app mesh...
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
CONSTANT
Real value(unsigned n, unsigned alpha, unsigned beta, Real x)
AuxVariableName _to_var_name
Variable in the main application to fill with the interpolation of the postprocessors.
SimpleRange< IndexType > as_range(const std::pair< IndexType, IndexType > &p)
MooseMesh wraps a libMesh::Mesh object and enhances its capabilities by caching additional data and s...
Definition: MooseMesh.h:88
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:31
void paramError(const std::string &param, Args... args) const
Emits an error prefixed with the file and line number of the given param (from the input file) along ...
static InputParameters validParams()
MONOMIAL
Real _power
Exponent for power-law decrease of the interpolation coefficients.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Base class for all MultiAppTransfer objects.
std::vector< std::unique_ptr< MultiAppCoordTransform > > _from_transforms
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
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...
unsigned int _num_points
Number of points to consider for the interpolation.
static System * find_sys(EquationSystems &es, const std::string &var_name)
Small helper function for finding the system containing the variable.
Definition: Transfer.C:89
PostprocessorName _postprocessor
Postprocessor in the child apps to transfer the data from.
MultiAppPostprocessorInterpolationTransfer(const InputParameters &parameters)
void ErrorVector unsigned int
bool _nodal
Whether the target variable is nodal or elemental.
uint8_t dof_id_type