www.mooseframework.org
MultiAppVariableValueSampleTransfer.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 "MooseVariable.h"
22 #include "MultiApp.h"
23 #include "SystemBase.h"
24 
25 #include "libmesh/meshfree_interpolation.h"
26 #include "libmesh/numeric_vector.h"
27 #include "libmesh/system.h"
28 
29 template <>
32 {
34  params.addRequiredParam<AuxVariableName>(
35  "variable", "The auxiliary variable to store the transferred values in.");
36  params.addRequiredParam<VariableName>("source_variable", "The variable to transfer from.");
37  return params;
38 }
39 
41  const InputParameters & parameters)
42  : MultiAppTransfer(parameters),
43  _to_var_name(getParam<AuxVariableName>("variable")),
44  _from_var_name(getParam<VariableName>("source_variable"))
45 {
46 }
47 
48 void
50 {
52 
53  _multi_app->problemBase().mesh().errorIfDistributedMesh("MultiAppVariableValueSampleTransfer");
54 }
55 
56 void
58 {
59  _console << "Beginning VariableValueSampleTransfer " << name() << std::endl;
60 
61  switch (_direction)
62  {
63  case TO_MULTIAPP:
64  {
65  FEProblemBase & from_problem = _multi_app->problemBase();
66  MooseVariable & from_var = from_problem.getVariable(0, _from_var_name);
67  SystemBase & from_system_base = from_var.sys();
68  SubProblem & from_sub_problem = from_system_base.subproblem();
69 
70  MooseMesh & from_mesh = from_problem.mesh();
71 
72  std::unique_ptr<PointLocatorBase> pl = from_mesh.getPointLocator();
73 
74  for (unsigned int i = 0; i < _multi_app->numGlobalApps(); i++)
75  {
76  Real value = -std::numeric_limits<Real>::max();
77 
78  { // Get the value of the variable at the point where this multiapp is in the master domain
79 
80  Point multi_app_position = _multi_app->position(i);
81 
82  std::vector<Point> point_vec(1, multi_app_position);
83 
84  // First find the element the hit lands in
85  const Elem * elem = (*pl)(multi_app_position);
86 
87  if (elem && elem->processor_id() == from_mesh.processor_id())
88  {
89  from_sub_problem.reinitElemPhys(elem, point_vec, 0);
90 
91  mooseAssert(from_var.sln().size() == 1, "No values in u!");
92  value = from_var.sln()[0];
93  }
94 
95  _communicator.max(value);
96 
97  if (value == -std::numeric_limits<Real>::max())
98  mooseError("Transfer failed to sample point value at point: ", multi_app_position);
99  }
100 
101  if (_multi_app->hasLocalApp(i))
102  {
103  Moose::ScopedCommSwapper swapper(_multi_app->comm());
104 
105  // Loop over the master nodes and set the value of the variable
106  System * to_sys = find_sys(_multi_app->appProblemBase(i).es(), _to_var_name);
107 
108  unsigned int sys_num = to_sys->number();
109  unsigned int var_num = to_sys->variable_number(_to_var_name);
110 
111  NumericVector<Real> & solution = _multi_app->appTransferVector(i, _to_var_name);
112 
113  MooseMesh & mesh = _multi_app->appProblemBase(i).mesh();
114 
115  MeshBase::const_node_iterator node_it = mesh.localNodesBegin();
116  MeshBase::const_node_iterator node_end = mesh.localNodesEnd();
117 
118  for (; node_it != node_end; ++node_it)
119  {
120  Node * node = *node_it;
121 
122  if (node->n_dofs(sys_num, var_num) > 0) // If this variable has dofs at this node
123  {
124  // The zero only works for LAGRANGE!
125  dof_id_type dof = node->dof_number(sys_num, var_num, 0);
126 
127  solution.set(dof, value);
128  }
129  }
130  solution.close();
131  _multi_app->appProblemBase(i).es().update();
132  }
133  }
134 
135  break;
136  }
137  case FROM_MULTIAPP:
138  {
139  mooseError("Doesn't make sense to transfer a sampled variable's value from a MultiApp!!");
140  break;
141  }
142  }
143 
144  _console << "Finished VariableValueSampleTransfer " << name() << std::endl;
145 }
const std::string & name() const
Get the name of the object.
Definition: MooseObject.h:47
void variableIntegrityCheck(const AuxVariableName &var_name) const
Utility to verify that the vEariable in the destination system exists.
Class for stuff related to variables.
Definition: MooseVariable.h:43
MeshBase::const_node_iterator localNodesBegin()
Calls local_nodes_begin/end() on the underlying libMesh mesh object.
Definition: MooseMesh.C:2000
MultiAppVariableValueSampleTransfer(const InputParameters &parameters)
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
Base class for a system (of equations)
Definition: SystemBase.h:91
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...
MeshBase::const_node_iterator localNodesEnd()
Definition: MooseMesh.C:2006
InputParameters validParams< MultiAppVariableValueSampleTransfer >()
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
virtual void initialSetup() override
Method called at the beginning of the simulation for checking integrity or doing one-time setup...
virtual SubProblem & subproblem()
Definition: SystemBase.h:103
virtual void reinitElemPhys(const Elem *elem, std::vector< Point > phys_points_in_elem, THREAD_ID tid)=0
virtual std::unique_ptr< PointLocatorBase > getPointLocator() const
Proxy function to get a (sub)PointLocator from either the underlying libMesh mesh (default)...
Definition: MooseMesh.C:2540
Generic class for solving transient nonlinear problems.
Definition: SubProblem.h:53
InputParameters validParams< MultiAppTransfer >()
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:250
Base class for all MultiAppTransfer objects.
const VariableValue & sln()
virtual MooseMesh & mesh() override
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
SystemBase & sys()
Get the system this variable is part of.
virtual MooseVariable & getVariable(THREAD_ID tid, const std::string &var_name) override
Returns the variable reference for requested variable which may be in any system. ...
virtual void execute() override
Execute the transfer.