www.mooseframework.org
DTKInterpolationHelper.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 #include "libmesh/libmesh_config.h"
15 
16 #ifdef LIBMESH_TRILINOS_HAVE_DTK
17 
18 // Moose Includes
19 #include "MooseError.h"
20 #include "DTKInterpolationHelper.h"
21 
22 #include "libmesh/equation_systems.h"
23 
24 // Ignore warnings coming from Trilinos/DTK headers.
25 #include "libmesh/ignore_warnings.h"
26 
27 // Trilinos Includes
28 #include <Teuchos_RCP.hpp>
29 #include <Teuchos_GlobalMPISession.hpp>
30 #include <Teuchos_Ptr.hpp>
31 #include <Teuchos_OpaqueWrapper.hpp>
32 
33 // DTK Includes
34 #include <DTK_MeshManager.hpp>
35 #include <DTK_MeshContainer.hpp>
36 #include <DTK_FieldManager.hpp>
37 #include <DTK_FieldContainer.hpp>
38 #include <DTK_FieldTools.hpp>
39 #include <DTK_CommTools.hpp>
40 #include <DTK_CommIndexer.hpp>
41 
42 // Restore warnings.
43 #include "libmesh/restore_warnings.h"
44 
45 namespace libMesh
46 {
47 
49 
51 {
52  for (std::map<EquationSystems *, DTKInterpolationAdapter *>::iterator it = adapters.begin();
53  it != adapters.end();
54  ++it)
55  delete it->second;
56 
57  for (std::map<std::pair<unsigned int, unsigned int>, shared_domain_map_type *>::iterator it =
58  dtk_maps.begin();
59  it != dtk_maps.end();
60  ++it)
61  delete it->second;
62 }
63 
64 void
66  unsigned int to,
67  const Variable * from_var,
68  const Variable * to_var,
69  const Point & from_offset,
70  const Point & to_offset,
71  MPI_Comm * from_mpi_comm,
72  MPI_Comm * to_mpi_comm)
73 {
74  Teuchos::RCP<const Teuchos::MpiComm<int>> from_comm;
75  Teuchos::RCP<const Teuchos::MpiComm<int>> to_comm;
76 
77  if (from_mpi_comm)
78  from_comm = Teuchos::rcp(new Teuchos::MpiComm<int>(
79  Teuchos::rcp(new Teuchos::OpaqueWrapper<MPI_Comm>(*from_mpi_comm))));
80 
81  if (to_mpi_comm)
82  to_comm = Teuchos::rcp(new Teuchos::MpiComm<int>(
83  Teuchos::rcp(new Teuchos::OpaqueWrapper<MPI_Comm>(*to_mpi_comm))));
84 
85  // Create a union comm for the shared domain transfer
86  Teuchos::RCP<const Teuchos::Comm<int>> comm_union;
87 
88  DataTransferKit::CommTools::unite(from_comm, to_comm, comm_union);
89 
90  EquationSystems * from_es = NULL;
91  EquationSystems * to_es = NULL;
92 
93  if (from_mpi_comm)
94  from_es = &from_var->system()->get_equation_systems();
95 
96  if (to_mpi_comm)
97  to_es = &to_var->system()->get_equation_systems();
98 
99  unsigned int dim = 3;
100 
101  if (from_es && to_es &&
102  (from_es->get_mesh().mesh_dimension() < to_es->get_mesh().mesh_dimension()))
103  mooseError(
104  "Receiving system dimension should be less than or equal to the sending system dimension!");
105 
106  if (from_es && to_es)
107  dim = std::max(from_es->get_mesh().mesh_dimension(), to_es->get_mesh().mesh_dimension());
108  else if (from_es)
109  dim = from_es->get_mesh().mesh_dimension();
110  else
111  dim = to_es->get_mesh().mesh_dimension();
112 
113  // Possibly make an Adapter for from_es
114  if (from_mpi_comm && adapters.find(from_es) == adapters.end())
115  adapters[from_es] = new DTKInterpolationAdapter(from_comm, *from_es, from_offset, dim);
116 
117  // Possibly make an Adapter for to_es
118  if (to_mpi_comm && adapters.find(to_es) == adapters.end())
119  adapters[to_es] = new DTKInterpolationAdapter(to_comm, *to_es, to_offset, dim);
120 
121  DTKInterpolationAdapter * from_adapter = NULL;
122  DTKInterpolationAdapter * to_adapter = NULL;
123 
124  if (from_mpi_comm)
125  from_adapter = adapters[from_es];
126 
127  if (to_mpi_comm)
128  to_adapter = adapters[to_es];
129 
130  std::pair<int, int> from_to(from, to);
131 
132  // If we haven't created a map for this pair of EquationSystems yet, do it now
133  if (dtk_maps.find(from_to) == dtk_maps.end())
134  {
135  shared_domain_map_type * map = NULL;
136 
137  if (from_mpi_comm)
138  map = new shared_domain_map_type(comm_union, from_es->get_mesh().mesh_dimension(), true);
139  else if (to_mpi_comm)
140  map = new shared_domain_map_type(comm_union, to_es->get_mesh().mesh_dimension(), true);
141 
142  dtk_maps[from_to] = map;
143 
144  // The tolerance here is for the "contains_point()" implementation in DTK. Set a larger value
145  // for a looser tolerance...
146  if (from_mpi_comm && to_mpi_comm)
147  {
148  if (to_var->type() == FEType())
149  map->setup(from_adapter->get_mesh_manager(),
150  to_adapter->get_target_coords(),
152  else
153  map->setup(from_adapter->get_mesh_manager(),
154  to_adapter->get_elem_target_coords(),
156  }
157  else if (from_mpi_comm)
158  map->setup(
159  from_adapter->get_mesh_manager(),
160  Teuchos::RCP<DataTransferKit::FieldManager<DTKInterpolationAdapter::MeshContainerType>>(),
162  else if (to_mpi_comm)
163  {
164  if (to_var->type() == FEType())
165  map->setup(Teuchos::RCP<
166  DataTransferKit::MeshManager<DTKInterpolationAdapter::MeshContainerType>>(),
167  to_adapter->get_target_coords(),
169  else
170  map->setup(Teuchos::RCP<
171  DataTransferKit::MeshManager<DTKInterpolationAdapter::MeshContainerType>>(),
172  to_adapter->get_elem_target_coords(),
174  }
175  }
176 
178 
179  if (from_mpi_comm)
180  from_evaluator = from_adapter->get_variable_evaluator(from_var->name());
181 
182  Teuchos::RCP<DataTransferKit::FieldManager<DTKInterpolationAdapter::FieldContainerType>>
183  to_values;
184 
185  if (to_mpi_comm)
186  to_values = to_adapter->get_values_to_fill(to_var->name());
187 
188  dtk_maps[from_to]->apply(from_evaluator, to_values);
189 
190  if (to_mpi_comm)
191  to_adapter->update_variable_values(to_var->name(), dtk_maps[from_to]->getMissedTargetPoints());
192 }
193 
194 } // namespace libMesh
195 
196 #endif // #ifdef LIBMESH_TRILINOS_HAVE_DTK
int eps(unsigned int i, unsigned int j)
2D version
void update_variable_values(std::string var_name, Teuchos::ArrayView< GlobalOrdinal > missed_points)
After computing values for a variable in this EquationSystems we need to take those values and put th...
RCP_Evaluator get_variable_evaluator(std::string var_name)
void mooseError(Args &&...args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:182
Teuchos::RCP< EvaluatorType > RCP_Evaluator
Teuchos::RCP< DataTransferKit::MeshManager< MeshContainerType > > get_mesh_manager()
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
DataTransferKit::SharedDomainMap< DTKInterpolationAdapter::MeshContainerType, DTKInterpolationAdapter::MeshContainerType > shared_domain_map_type
void transferWithOffset(unsigned int from, unsigned int to, const Variable *from_var, const Variable *to_var, const Point &from_offset, const Point &to_offset, MPI_Comm *from_mpi_comm, MPI_Comm *to_mpi_comm)
Do an interpolation with possibly offsetting each of the domains (moving them).
std::map< EquationSystems *, DTKInterpolationAdapter * > adapters
The DTKAdapter associated with each EquationSystems.
Teuchos::RCP< DataTransferKit::FieldManager< MeshContainerType > > get_elem_target_coords()
Used to get the centroids for the receiving elements.
Teuchos::RCP< DataTransferKit::FieldManager< FieldContainerType > > get_values_to_fill(std::string var_name)
Teuchos::RCP< DataTransferKit::FieldManager< MeshContainerType > > get_target_coords()
std::map< std::pair< unsigned int, unsigned int >, shared_domain_map_type * > dtk_maps
The dtk shared domain maps for pairs of EquationSystems (from, to)