20 #include "libmesh/libmesh_config.h" 22 #ifdef LIBMESH_TRILINOS_HAVE_DTK 25 #include "libmesh/dtk_adapter.h" 26 #include "libmesh/dtk_evaluator.h" 27 #include "libmesh/int_range.h" 28 #include "libmesh/mesh.h" 29 #include "libmesh/numeric_vector.h" 30 #include "libmesh/elem.h" 31 #include "libmesh/equation_systems.h" 34 #include "libmesh/ignore_warnings.h" 35 #include <DTK_MeshTypes.hpp> 36 #include <Teuchos_Comm.hpp> 37 #include "libmesh/restore_warnings.h" 49 mesh(in_es.get_mesh()),
52 std::set<unsigned int> semi_local_nodes;
64 for (
const auto &
id : semi_local_nodes)
70 for (
unsigned int j=0; j<
dim; j++)
83 Teuchos::ArrayRCP<int> elements(n_local_elem);
84 Teuchos::ArrayRCP<int> connectivity(n_nodes_per_elem*n_local_elem);
90 for (
const auto & elem :
as_range(
mesh.local_elements_begin(),
91 mesh.local_elements_end()))
93 elements[i] = elem->id();
95 for (
unsigned int j=0; j<n_nodes_per_elem; j++)
96 connectivity[(j*n_local_elem)+i] = elem->node_id(j);
102 Teuchos::ArrayRCP<int> permutation_list(n_nodes_per_elem);
103 std::iota(permutation_list.begin(), permutation_list.end(), 0);
154 Teuchos::RCP<MeshContainerType>
165 Teuchos::ArrayRCP<Teuchos::RCP<MeshContainerType>> mesh_blocks(1);
166 mesh_blocks[0] = mesh_container;
169 mesh_manager = Teuchos::rcp(
new DataTransferKit::MeshManager<MeshContainerType>(mesh_blocks,
comm,
dim) );
172 target_coords = Teuchos::rcp(
new DataTransferKit::FieldManager<MeshContainerType>(mesh_container,
comm));
179 auto [it, emplaced] =
evaluators.emplace(var_name,
nullptr);
186 it->second = Teuchos::rcp(
new DTKEvaluator(*sys, var_name));
192 Teuchos::RCP<DataTransferKit::FieldManager<DTKAdapter::FieldContainerType>>
203 Teuchos::RCP<FieldContainerType> field_container = Teuchos::rcp(
new FieldContainerType(data_space, 1));
204 it->second = Teuchos::rcp(
new DataTransferKit::FieldManager<FieldContainerType>(field_container,
comm));
216 Teuchos::RCP<FieldContainerType> values =
values_to_fill[var_name]->field();
220 for (
const auto &
value : *values)
222 unsigned int node_num =
vertices[i];
264 DataTransferKit::DTK_ElementTopology
270 return DataTransferKit::DTK_LINE_SEGMENT;
271 else if (type ==
TRI3)
272 return DataTransferKit::DTK_TRIANGLE;
273 else if (type ==
QUAD4)
274 return DataTransferKit::DTK_QUADRILATERAL;
275 else if (type ==
TET4)
276 return DataTransferKit::DTK_TETRAHEDRON;
277 else if (type ==
HEX8)
278 return DataTransferKit::DTK_HEXAHEDRON;
280 return DataTransferKit::DTK_PYRAMID;
282 libmesh_error_msg(
"Element type not supported by DTK!");
288 for (
const auto & elem :
as_range(
mesh.local_elements_begin(),
mesh.local_elements_end()))
289 for (
const Node & node : elem->node_ref_range())
290 semi_local_nodes.insert(node.id());
295 #endif // #ifdef LIBMESH_TRILINOS_HAVE_DTK ElemType
Defines an enum for geometric element types.
This is the EquationSystems class.
dof_id_type dof_number(const unsigned int s, const unsigned int var, const unsigned int comp) const
std::map< std::string, RCP_Evaluator > evaluators
Map of variable names to RCP_Evaluator objects.
A Node is like a Point, but with more information.
unsigned int n_systems() const
Teuchos::RCP< const Teuchos::Comm< int > > comm
void update_variable_values(std::string var_name)
After computing values for a variable in this EquationSystems we need to take those values and put th...
unsigned int num_local_nodes
This is the base class from which all geometric element types are derived.
The libMesh namespace provides an interface to certain functionality in the library.
Implements the evaluate() function to compute FE solution values at points requested by DTK...
const T_sys & get_system(std::string_view name) const
DataTransferKit::MeshContainer< int > MeshContainerType
dof_id_type n_local_elem() const
void iota(ForwardIter first, ForwardIter last, T value)
Utility::iota was created back when std::iota was just an SGI STL extension.
DataTransferKit::DTK_ElementTopology get_element_topology(const Elem *elem)
unsigned int variable_number(std::string_view var) const
unsigned int number() const
RCP_Evaluator get_variable_evaluator(std::string var_name)
virtual unsigned int n_nodes() const =0
Manages consistently variables, degrees of freedom, and coefficient vectors.
SimpleRange< IndexType > as_range(const std::pair< IndexType, IndexType > &p)
Helper function that allows us to treat a homogenous pair as a range.
std::unique_ptr< NumericVector< Number > > solution
Data structure to hold solution values.
DTKAdapter(Teuchos::RCP< const Teuchos::Comm< int >> in_comm, EquationSystems &in_es)
Teuchos::RCP< EvaluatorType > RCP_Evaluator
virtual const Elem * elem_ptr(const dof_id_type i) const =0
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
std::map< std::string, Teuchos::RCP< DataTransferKit::FieldManager< FieldContainerType > > > values_to_fill
Map of variable names to arrays to be filled by a transfer.
Teuchos::RCP< DataTransferKit::FieldManager< MeshContainerType > > target_coords
void get_semi_local_nodes(std::set< unsigned int > &semi_local_nodes)
Helper function that fills the std::set with all of the node numbers of nodes connected to local elem...
virtual const Node & node_ref(const dof_id_type i) const
processor_id_type processor_id() const
Teuchos::RCP< DataTransferKit::MeshManager< MeshContainerType > > mesh_manager
Teuchos::ArrayRCP< int > vertices
processor_id_type processor_id() const
virtual ElemType type() const =0
DataTransferKit::FieldContainer< double > FieldContainerType
System * find_sys(std::string var_name)
Small helper function for finding the system containing the variable.
Teuchos::RCP< DataTransferKit::FieldManager< FieldContainerType > > get_values_to_fill(std::string var_name)