19 #include "libmesh/meshfree_interpolation.h" 20 #include "libmesh/numeric_vector.h" 21 #include "libmesh/system.h" 22 #include "libmesh/radial_basis_interpolation.h" 30 params.
addClassDescription(
"Transfer postprocessor data from sub-application into field data on " 31 "the parent application.");
33 "variable",
"The auxiliary variable to store the transferred values in.");
34 params.
addRequiredParam<PostprocessorName>(
"postprocessor",
"The Postprocessor to interpolate.");
36 "num_points", 3,
"The number of nearest points to use for interpolation.");
38 "power", 2,
"The polynomial power to use for calculation of the decay in the interpolation.");
40 MooseEnum interp_type(
"inverse_distance radial_basis",
"inverse_distance");
42 params.
addParam<
MooseEnum>(
"interp_type", interp_type,
"The algorithm to use for interpolation.");
46 "Radius to use for radial_basis interpolation. If negative " 47 "then the radius is taken as the max distance between " 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")),
64 paramError(
"to_multi_app",
"Unused parameter; only from-MultiApp transfers are implemented");
70 paramError(
"variable",
"Must be either CONSTANT MONOMIAL or FIRST LAGRANGE");
78 TIME_SECTION(
"MultiAppPostprocessorInterpolationTransfer::execute()",
80 "Transferring/interpolating postprocessors");
86 mooseError(
"Interpolation from a variable to a MultiApp postprocessors has not been " 87 "implemented. Use MultiAppVariableValueSamplePostprocessorTransfer!");
92 std::unique_ptr<InverseDistanceInterpolation<LIBMESH_DIM>> idi;
97 idi = std::make_unique<InverseDistanceInterpolation<LIBMESH_DIM>>(
107 std::vector<Point> & src_pts(idi->get_source_points());
108 std::vector<Number> & src_vals(idi->get_source_vals());
110 std::vector<std::string> field_vars;
112 idi->set_field_variables(field_vars);
115 mooseAssert(
_to_transforms.size() == 1,
"There should only be one transform here");
124 src_pts.push_back(to_coord_transform.mapBack((*
_from_transforms[i])(Point(0))));
128 src_pts.push_back(to_coord_transform.mapBack(
136 idi->prepare_for_use();
142 unsigned int sys_num = to_sys->number();
143 unsigned int var_num = to_sys->variable_number(
_to_var_name);
145 NumericVector<Real> & solution = *to_sys->solution;
149 std::vector<std::string> vars;
156 for (
const auto & node :
as_range(
mesh.localNodesBegin(),
mesh.localNodesEnd()))
158 if (node->n_dofs(sys_num, var_num) > 0)
160 std::vector<Point> pts;
161 std::vector<Number> vals;
163 pts.push_back(*node);
166 idi->interpolate_field_data(vars, pts, vals);
171 dof_id_type dof = node->dof_number(sys_num, var_num, 0);
173 solution.set(dof,
value);
180 for (
const auto & elem :
181 as_range(
mesh.getMesh().local_elements_begin(),
mesh.getMesh().local_elements_end()))
184 if (elem->n_dofs(sys_num, var_num) > 0)
186 std::vector<Point> pts;
187 std::vector<Number> vals;
189 pts.push_back(elem->vertex_average());
192 idi->interpolate_field_data(vars, pts, vals);
196 dof_id_type dof = elem->dof_number(sys_num, var_num, 0);
197 solution.set(dof,
value);
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
registerMooseObject("MooseApp", MultiAppPostprocessorInterpolationTransfer)
Real _radius
Radius to consider for inverse interpolation.
const Parallel::Communicator & _communicator
MooseEnum _interp_type
Type of interpolation method.
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.
Real value(unsigned n, unsigned alpha, unsigned beta, Real x)
virtual void execute() override
Execute the transfer.
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...
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
void paramError(const std::string ¶m, 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()
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.
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.
PostprocessorName _postprocessor
Postprocessor in the child apps to transfer the data from.
static InputParameters validParams()
MultiAppPostprocessorInterpolationTransfer(const InputParameters ¶meters)
void ErrorVector unsigned int
bool _nodal
Whether the target variable is nodal or elemental.