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