www.mooseframework.org
DTKInterpolationAdapter.h
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #pragma once
11 
12 #include "libmesh/libmesh_config.h"
13 
14 #ifdef LIBMESH_TRILINOS_HAVE_DTK
15 
16 #include "libmesh/point.h"
17 
18 // Ignore warnings coming from Trilinos/DTK.
19 #include "libmesh/ignore_warnings.h"
20 
21 // DTK includes
22 #include <DTK_MeshManager.hpp>
23 #include <DTK_MeshContainer.hpp>
24 #include <DTK_MeshTraits.hpp>
25 #include <DTK_MeshTraitsFieldAdapter.hpp>
26 #include <DTK_FieldManager.hpp>
27 #include <DTK_FieldContainer.hpp>
28 #include <DTK_FieldEvaluator.hpp>
29 
30 // Trilinos includes
31 #include <Teuchos_RCP.hpp>
32 #include <Teuchos_ArrayRCP.hpp>
33 #include <Teuchos_DefaultMpiComm.hpp>
34 
35 // Restore warnings.
36 #include "libmesh/restore_warnings.h"
37 
38 // Forward declarations
39 namespace libMesh
40 {
41 class EquationSystems;
42 class Elem;
43 class MeshBase;
44 class System;
45 }
46 
48 {
49 public:
50  DTKInterpolationAdapter(Teuchos::RCP<const Teuchos::MpiComm<int>> in_comm,
51  libMesh::EquationSystems & in_es,
52  const libMesh::Point & offset,
53  unsigned int from_dim);
54 
55  typedef DataTransferKit::MeshContainer<long unsigned int> MeshContainerType;
56  typedef DataTransferKit::FieldContainer<double> FieldContainerType;
57  typedef DataTransferKit::MeshTraits<MeshContainerType>::global_ordinal_type GlobalOrdinal;
58  typedef DataTransferKit::FieldEvaluator<GlobalOrdinal, FieldContainerType> EvaluatorType;
59  typedef Teuchos::RCP<EvaluatorType> RCP_Evaluator;
60 
61  Teuchos::RCP<DataTransferKit::MeshManager<MeshContainerType>> get_mesh_manager()
62  {
63  return mesh_manager;
64  }
65  RCP_Evaluator get_variable_evaluator(std::string var_name);
66  Teuchos::RCP<DataTransferKit::FieldManager<MeshContainerType>> get_target_coords()
67  {
68  return target_coords;
69  }
70 
74  Teuchos::RCP<DataTransferKit::FieldManager<MeshContainerType>> get_elem_target_coords()
75  {
76  return elem_centroid_coords;
77  }
78 
79  Teuchos::RCP<DataTransferKit::FieldManager<FieldContainerType>>
80  get_values_to_fill(std::string var_name);
81 
86  void update_variable_values(std::string var_name,
87  Teuchos::ArrayView<GlobalOrdinal> missed_points);
88 
89 protected:
93  DataTransferKit::DTK_ElementTopology get_element_topology(const libMesh::Elem * elem);
94 
99  void get_semi_local_nodes(std::set<GlobalOrdinal> & semi_local_nodes);
100 
101  Teuchos::RCP<const Teuchos::MpiComm<int>> comm;
102  libMesh::EquationSystems & es;
103  libMesh::Point _offset;
104  const libMesh::MeshBase & mesh;
105  unsigned int dim;
106 
107  unsigned int num_local_nodes;
108  Teuchos::ArrayRCP<GlobalOrdinal> vertices;
109  Teuchos::ArrayRCP<GlobalOrdinal> elements;
110 
111  Teuchos::RCP<DataTransferKit::MeshManager<MeshContainerType>> mesh_manager;
113  Teuchos::RCP<DataTransferKit::FieldManager<MeshContainerType>> target_coords;
114  Teuchos::RCP<DataTransferKit::FieldManager<MeshContainerType>> elem_centroid_coords;
115 
117  std::map<std::string, Teuchos::RCP<DataTransferKit::FieldManager<FieldContainerType>>>
119 
121  std::map<std::string, RCP_Evaluator> evaluators;
122 };
123 
124 #endif // #ifdef LIBMESH_TRILINOS_HAVE_DTK
125 
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)
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
void get_semi_local_nodes(std::set< GlobalOrdinal > &semi_local_nodes)
Helper function that fills the std::set with all of the node numbers of nodes connected to local elem...
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::DTK_ElementTopology get_element_topology(const libMesh::Elem *elem)
Helper that returns the DTK ElementTopology for a given Elem.
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< FieldContainerType > > get_values_to_fill(std::string var_name)
Teuchos::RCP< DataTransferKit::FieldManager< MeshContainerType > > get_target_coords()
DTKInterpolationAdapter(Teuchos::RCP< const Teuchos::MpiComm< int >> in_comm, libMesh::EquationSystems &in_es, const libMesh::Point &offset, unsigned int from_dim)