www.mooseframework.org
DTKInterpolationAdapter.h
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 
15 #ifndef DTKINTERPOLATIONADAPTER_H
16 #define DTKINTERPOLATIONADAPTER_H
17 
18 #include "libmesh/libmesh_config.h"
19 
20 #ifdef LIBMESH_TRILINOS_HAVE_DTK
21 
22 #include "libmesh/point.h"
23 
24 // Ignore warnings coming from Trilinos/DTK.
25 #include "libmesh/ignore_warnings.h"
26 
27 // DTK includes
28 #include <DTK_MeshManager.hpp>
29 #include <DTK_MeshContainer.hpp>
30 #include <DTK_MeshTraits.hpp>
31 #include <DTK_MeshTraitsFieldAdapter.hpp>
32 #include <DTK_FieldManager.hpp>
33 #include <DTK_FieldContainer.hpp>
34 #include <DTK_FieldEvaluator.hpp>
35 
36 // Trilinos includes
37 #include <Teuchos_RCP.hpp>
38 #include <Teuchos_ArrayRCP.hpp>
39 #include <Teuchos_DefaultMpiComm.hpp>
40 
41 // Restore warnings.
42 #include "libmesh/restore_warnings.h"
43 
44 // Forward declarations
45 namespace libMesh
46 {
47 class EquationSystems;
48 class Elem;
49 class MeshBase;
50 class System;
51 }
52 
54 {
55 public:
56  DTKInterpolationAdapter(Teuchos::RCP<const Teuchos::MpiComm<int>> in_comm,
57  libMesh::EquationSystems & in_es,
58  const libMesh::Point & offset,
59  unsigned int from_dim);
60 
61  typedef DataTransferKit::MeshContainer<long unsigned int> MeshContainerType;
62  typedef DataTransferKit::FieldContainer<double> FieldContainerType;
63  typedef DataTransferKit::MeshTraits<MeshContainerType>::global_ordinal_type GlobalOrdinal;
64  typedef DataTransferKit::FieldEvaluator<GlobalOrdinal, FieldContainerType> EvaluatorType;
65  typedef Teuchos::RCP<EvaluatorType> RCP_Evaluator;
66 
67  Teuchos::RCP<DataTransferKit::MeshManager<MeshContainerType>> get_mesh_manager()
68  {
69  return mesh_manager;
70  }
71  RCP_Evaluator get_variable_evaluator(std::string var_name);
72  Teuchos::RCP<DataTransferKit::FieldManager<MeshContainerType>> get_target_coords()
73  {
74  return target_coords;
75  }
76 
80  Teuchos::RCP<DataTransferKit::FieldManager<MeshContainerType>> get_elem_target_coords()
81  {
82  return elem_centroid_coords;
83  }
84 
85  Teuchos::RCP<DataTransferKit::FieldManager<FieldContainerType>>
86  get_values_to_fill(std::string var_name);
87 
92  void update_variable_values(std::string var_name,
93  Teuchos::ArrayView<GlobalOrdinal> missed_points);
94 
95 protected:
99  DataTransferKit::DTK_ElementTopology get_element_topology(const libMesh::Elem * elem);
100 
105  void get_semi_local_nodes(std::set<GlobalOrdinal> & semi_local_nodes);
106 
107  Teuchos::RCP<const Teuchos::MpiComm<int>> comm;
108  libMesh::EquationSystems & es;
109  libMesh::Point _offset;
110  const libMesh::MeshBase & mesh;
111  unsigned int dim;
112 
113  unsigned int num_local_nodes;
114  Teuchos::ArrayRCP<GlobalOrdinal> vertices;
115  Teuchos::ArrayRCP<GlobalOrdinal> elements;
116 
117  Teuchos::RCP<DataTransferKit::MeshManager<MeshContainerType>> mesh_manager;
118  RCP_Evaluator field_evaluator;
119  Teuchos::RCP<DataTransferKit::FieldManager<MeshContainerType>> target_coords;
120  Teuchos::RCP<DataTransferKit::FieldManager<MeshContainerType>> elem_centroid_coords;
121 
123  std::map<std::string, Teuchos::RCP<DataTransferKit::FieldManager<FieldContainerType>>>
125 
127  std::map<std::string, RCP_Evaluator> evaluators;
128 };
129 
130 #endif // #ifdef LIBMESH_TRILINOS_HAVE_DTK
131 
132 #endif // #define DTKINTERPOLATIONADAPTER_H
std::map< std::string, Teuchos::RCP< DataTransferKit::FieldManager< FieldContainerType > > > values_to_fill
Map of variable names to arrays to be filled by a transfer.
DataTransferKit::FieldContainer< double > FieldContainerType
Teuchos::RCP< EvaluatorType > RCP_Evaluator
Teuchos::RCP< DataTransferKit::MeshManager< MeshContainerType > > get_mesh_manager()
DataTransferKit::FieldEvaluator< GlobalOrdinal, FieldContainerType > EvaluatorType
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
std::map< std::string, RCP_Evaluator > evaluators
Map of variable names to RCP_Evaluator objects.
DataTransferKit::MeshTraits< MeshContainerType >::global_ordinal_type GlobalOrdinal
Teuchos::ArrayRCP< GlobalOrdinal > vertices
Teuchos::RCP< DataTransferKit::MeshManager< MeshContainerType > > mesh_manager
DataTransferKit::MeshContainer< long unsigned int > MeshContainerType
const libMesh::MeshBase & mesh
Teuchos::RCP< DataTransferKit::FieldManager< MeshContainerType > > elem_centroid_coords
Teuchos::RCP< DataTransferKit::FieldManager< MeshContainerType > > target_coords
Teuchos::RCP< const Teuchos::MpiComm< int > > comm
libMesh::EquationSystems & es
Teuchos::RCP< DataTransferKit::FieldManager< MeshContainerType > > get_elem_target_coords()
Used to get the centroids for the receiving elements.
Teuchos::ArrayRCP< GlobalOrdinal > elements
Teuchos::RCP< DataTransferKit::FieldManager< MeshContainerType > > get_target_coords()