libMesh
Public Types | Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
libMesh::DTKAdapter Class Reference

The DTKAdapter is used with the DTKSolutionTransfer object to adapt libmesh data to the DTK interface. More...

#include <dtk_adapter.h>

Public Types

typedef DataTransferKit::MeshContainer< intMeshContainerType
 
typedef DataTransferKit::FieldContainer< double > FieldContainerType
 
typedef DataTransferKit::MeshTraits< MeshContainerType >::global_ordinal_type GlobalOrdinal
 
typedef DataTransferKit::FieldEvaluator< GlobalOrdinal, FieldContainerTypeEvaluatorType
 
typedef Teuchos::RCP< EvaluatorTypeRCP_Evaluator
 

Public Member Functions

 DTKAdapter (Teuchos::RCP< const Teuchos::Comm< int >> in_comm, EquationSystems &in_es)
 
Teuchos::RCP< DataTransferKit::MeshManager< MeshContainerType > > get_mesh_manager ()
 
RCP_Evaluator get_variable_evaluator (std::string var_name)
 
Teuchos::RCP< DataTransferKit::FieldManager< MeshContainerType > > get_target_coords ()
 
Teuchos::RCP< DataTransferKit::FieldManager< FieldContainerType > > get_values_to_fill (std::string var_name)
 
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 them in the actual solution vector. More...
 

Protected Member Functions

Systemfind_sys (std::string var_name)
 Small helper function for finding the system containing the variable. More...
 
DataTransferKit::DTK_ElementTopology get_element_topology (const Elem *elem)
 
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 elements. More...
 

Protected Attributes

Teuchos::RCP< const Teuchos::Comm< int > > comm
 
EquationSystemses
 
const MeshBasemesh
 
unsigned int dim
 
unsigned int num_local_nodes
 
Teuchos::ArrayRCP< intvertices
 
Teuchos::RCP< DataTransferKit::MeshManager< MeshContainerType > > mesh_manager
 
RCP_Evaluator field_evaluator
 
Teuchos::RCP< DataTransferKit::FieldManager< MeshContainerType > > target_coords
 
std::map< std::string, Teuchos::RCP< DataTransferKit::FieldManager< FieldContainerType > > > values_to_fill
 Map of variable names to arrays to be filled by a transfer. More...
 
std::map< std::string, RCP_Evaluatorevaluators
 Map of variable names to RCP_Evaluator objects. More...
 

Detailed Description

The DTKAdapter is used with the DTKSolutionTransfer object to adapt libmesh data to the DTK interface.

Author
Derek Gaston
Date
2013

Definition at line 52 of file dtk_adapter.h.

Member Typedef Documentation

typedef DataTransferKit::FieldEvaluator<GlobalOrdinal,FieldContainerType> libMesh::DTKAdapter::EvaluatorType

Definition at line 60 of file dtk_adapter.h.

typedef DataTransferKit::FieldContainer<double> libMesh::DTKAdapter::FieldContainerType

Definition at line 58 of file dtk_adapter.h.

typedef DataTransferKit::MeshTraits<MeshContainerType>::global_ordinal_type libMesh::DTKAdapter::GlobalOrdinal

Definition at line 59 of file dtk_adapter.h.

typedef DataTransferKit::MeshContainer<int> libMesh::DTKAdapter::MeshContainerType

Definition at line 57 of file dtk_adapter.h.

Definition at line 61 of file dtk_adapter.h.

Constructor & Destructor Documentation

libMesh::DTKAdapter::DTKAdapter ( Teuchos::RCP< const Teuchos::Comm< int >>  in_comm,
EquationSystems in_es 
)

Definition at line 42 of file dtk_adapter.C.

References comm, dim, libMesh::MeshBase::elem(), end, get_element_topology(), get_semi_local_nodes(), libMesh::DofObject::id(), libMesh::MeshBase::local_elements_begin(), libMesh::MeshBase::local_elements_end(), mesh, mesh_manager, libMesh::MeshBase::n_local_elem(), libMesh::Elem::n_nodes(), libMesh::Elem::node_id(), libMesh::MeshBase::node_ref(), num_local_nodes, target_coords, and vertices.

42  :
43  comm(in_comm),
44  es(in_es),
45  mesh(in_es.get_mesh()),
47 {
48  std::set<unsigned int> semi_local_nodes;
49  get_semi_local_nodes(semi_local_nodes);
50 
51  num_local_nodes = semi_local_nodes.size();
52 
53  vertices.resize(num_local_nodes);
54  Teuchos::ArrayRCP<double> coordinates(num_local_nodes * dim);
55 
56  // Fill in the vertices and coordinates
57  {
58  unsigned int i = 0;
59 
60  for (std::set<unsigned int>::iterator it = semi_local_nodes.begin();
61  it != semi_local_nodes.end();
62  ++it)
63  {
64  const Node & node = mesh.node_ref(*it);
65 
66  vertices[i] = node.id();
67 
68  for (unsigned int j=0; j<dim; j++)
69  coordinates[(j*num_local_nodes) + i] = node(j);
70 
71  i++;
72  }
73  }
74 
75  // Currently assuming all elements are the same!
76  DataTransferKit::DTK_ElementTopology element_topology = get_element_topology(mesh.elem(0));
77  unsigned int n_nodes_per_elem = mesh.elem(0)->n_nodes();
78 
79  unsigned int n_local_elem = mesh.n_local_elem();
80 
81  Teuchos::ArrayRCP<int> elements(n_local_elem);
82  Teuchos::ArrayRCP<int> connectivity(n_nodes_per_elem*n_local_elem);
83 
84  // Fill in the elements and connectivity
85  {
86  unsigned int i = 0;
87 
88  MeshBase::const_element_iterator end = mesh.local_elements_end();
89  for (MeshBase::const_element_iterator it = mesh.local_elements_begin();
90  it != end;
91  ++it)
92  {
93  const Elem & elem = *(*it);
94  elements[i] = elem.id();
95 
96  for (unsigned int j=0; j<n_nodes_per_elem; j++)
97  connectivity[(j*n_local_elem)+i] = elem.node_id(j);
98 
99  i++;
100  }
101  }
102 
103  Teuchos::ArrayRCP<int> permutation_list(n_nodes_per_elem);
104  for (unsigned int i = 0; i < n_nodes_per_elem; ++i )
105  permutation_list[i] = i;
106 
107  /*
108  if (this->processor_id() == 1)
109  sleep(1);
110 
111  libMesh::out<<"n_nodes_per_elem: "<<n_nodes_per_elem<<std::endl;
112 
113  libMesh::out<<"Dim: "<<dim<<std::endl;
114 
115  libMesh::err<<"Vertices size: "<<vertices.size()<<std::endl;
116  {
117  libMesh::err<<this->processor_id()<<" Vertices: ";
118 
119  for (std::size_t i=0; i<vertices.size(); i++)
120  libMesh::err<<vertices[i]<<" ";
121 
122  libMesh::err<<std::endl;
123  }
124 
125  libMesh::err<<"Coordinates size: "<<coordinates.size()<<std::endl;
126  {
127  libMesh::err<<this->processor_id()<<" Coordinates: ";
128 
129  for (std::size_t i=0; i<coordinates.size(); i++)
130  libMesh::err<<coordinates[i]<<" ";
131 
132  libMesh::err<<std::endl;
133  }
134 
135  libMesh::err<<"Connectivity size: "<<connectivity.size()<<std::endl;
136  {
137  libMesh::err<<this->processor_id()<<" Connectivity: ";
138 
139  for (std::size_t i=0; i<connectivity.size(); i++)
140  libMesh::err<<connectivity[i]<<" ";
141 
142  libMesh::err<<std::endl;
143  }
144 
145  libMesh::err<<"Permutation_List size: "<<permutation_list.size()<<std::endl;
146  {
147  libMesh::err<<this->processor_id()<<" Permutation_List: ";
148 
149  for (std::size_t i=0; i<permutation_list.size(); i++)
150  libMesh::err<<permutation_list[i]<<" ";
151 
152  libMesh::err<<std::endl;
153  }
154 
155  */
156  Teuchos::RCP<MeshContainerType>
157  mesh_container = Teuchos::rcp(new MeshContainerType(dim,
158  vertices,
159  coordinates,
160  element_topology,
161  n_nodes_per_elem,
162  elements,
163  connectivity,
164  permutation_list));
165 
166  // We only have 1 element topology in this grid so we make just one mesh block
167  Teuchos::ArrayRCP<Teuchos::RCP<MeshContainerType>> mesh_blocks(1);
168  mesh_blocks[0] = mesh_container;
169 
170  // Create the MeshManager
171  mesh_manager = Teuchos::rcp(new DataTransferKit::MeshManager<MeshContainerType>(mesh_blocks, comm, dim) );
172 
173  // Pack the coordinates into a field, this will be the positions we'll ask for other systems fields at
174  target_coords = Teuchos::rcp(new DataTransferKit::FieldManager<MeshContainerType>(mesh_container, comm));
175 }
Teuchos::RCP< const Teuchos::Comm< int > > comm
Definition: dtk_adapter.h:94
unsigned int dim
Definition: dtk_adapter.h:97
virtual element_iterator local_elements_begin()=0
virtual const Node & node_ref(const dof_id_type i) const
Definition: mesh_base.h:420
unsigned int num_local_nodes
Definition: dtk_adapter.h:99
virtual const Elem * elem(const dof_id_type i) const
Definition: mesh_base.h:523
IterBase * end
Also have a polymorphic pointer to the end object, this prevents iterating past the end...
DataTransferKit::MeshContainer< int > MeshContainerType
Definition: dtk_adapter.h:57
const MeshBase & mesh
Definition: dtk_adapter.h:96
DataTransferKit::DTK_ElementTopology get_element_topology(const Elem *elem)
Definition: dtk_adapter.C:259
virtual unsigned int n_nodes() const =0
EquationSystems & es
Definition: dtk_adapter.h:95
virtual element_iterator local_elements_end()=0
dof_id_type n_local_elem() const
Definition: mesh_base.h:372
unsigned int mesh_dimension() const
Definition: mesh_base.C:148
Teuchos::RCP< DataTransferKit::FieldManager< MeshContainerType > > target_coords
Definition: dtk_adapter.h:104
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...
Definition: dtk_adapter.C:280
Teuchos::RCP< DataTransferKit::MeshManager< MeshContainerType > > mesh_manager
Definition: dtk_adapter.h:102
Teuchos::ArrayRCP< int > vertices
Definition: dtk_adapter.h:100

Member Function Documentation

System * libMesh::DTKAdapter::find_sys ( std::string  var_name)
protected

Small helper function for finding the system containing the variable.

Note
This implies that variable names are unique across all systems!

Note that this implies that variable names are unique across all systems!

Definition at line 239 of file dtk_adapter.C.

References es, libMesh::EquationSystems::get_system(), libMesh::libmesh_assert(), libmesh_nullptr, libMesh::EquationSystems::n_systems(), and libMesh::sys.

Referenced by get_target_coords(), get_variable_evaluator(), and update_variable_values().

240 {
241  System * sys = libmesh_nullptr;
242 
243  // Find the system this variable is from
244  for (unsigned int i=0; i<es.n_systems(); i++)
245  {
246  if (es.get_system(i).has_variable(var_name))
247  {
248  sys = &es.get_system(i);
249  break;
250  }
251  }
252 
253  libmesh_assert(sys);
254 
255  return sys;
256 }
ImplicitSystem & sys
const class libmesh_nullptr_t libmesh_nullptr
libmesh_assert(j)
EquationSystems & es
Definition: dtk_adapter.h:95
unsigned int n_systems() const
const T_sys & get_system(const std::string &name) const
DataTransferKit::DTK_ElementTopology libMesh::DTKAdapter::get_element_topology ( const Elem elem)
protected
Returns
The DTK ElementTopology for a given Elem.

Definition at line 259 of file dtk_adapter.C.

References libMesh::EDGE2, libMesh::HEX8, libMesh::PYRAMID5, libMesh::QUAD4, libMesh::TET4, libMesh::TRI3, and libMesh::Elem::type().

Referenced by DTKAdapter(), and get_target_coords().

260 {
261  ElemType type = elem->type();
262 
263  if (type == EDGE2)
264  return DataTransferKit::DTK_LINE_SEGMENT;
265  else if (type == TRI3)
266  return DataTransferKit::DTK_TRIANGLE;
267  else if (type == QUAD4)
268  return DataTransferKit::DTK_QUADRILATERAL;
269  else if (type == TET4)
270  return DataTransferKit::DTK_TETRAHEDRON;
271  else if (type == HEX8)
272  return DataTransferKit::DTK_HEXAHEDRON;
273  else if (type == PYRAMID5)
274  return DataTransferKit::DTK_PYRAMID;
275 
276  libmesh_error_msg("Element type not supported by DTK!");
277 }
ElemType
Defines an enum for geometric element types.
Teuchos::RCP<DataTransferKit::MeshManager<MeshContainerType> > libMesh::DTKAdapter::get_mesh_manager ( )

Definition at line 64 of file dtk_adapter.h.

References get_variable_evaluator(), and mesh_manager.

Referenced by libMesh::DTKSolutionTransfer::transfer().

64 { return mesh_manager; }
Teuchos::RCP< DataTransferKit::MeshManager< MeshContainerType > > mesh_manager
Definition: dtk_adapter.h:102
void libMesh::DTKAdapter::get_semi_local_nodes ( std::set< unsigned int > &  semi_local_nodes)
protected

Helper function that fills the std::set with all of the node numbers of nodes connected to local elements.

Definition at line 280 of file dtk_adapter.C.

References end, libMesh::MeshBase::local_elements_begin(), libMesh::MeshBase::local_elements_end(), mesh, libMesh::Elem::n_nodes(), and libMesh::Elem::node_id().

Referenced by DTKAdapter(), and get_target_coords().

281 {
282  MeshBase::const_element_iterator end = mesh.local_elements_end();
283  for (MeshBase::const_element_iterator it = mesh.local_elements_begin();
284  it != end;
285  ++it)
286  {
287  const Elem & elem = *(*it);
288 
289  for (unsigned int j=0; j<elem.n_nodes(); j++)
290  semi_local_nodes.insert(elem.node_id(j));
291  }
292 }
virtual element_iterator local_elements_begin()=0
IterBase * end
Also have a polymorphic pointer to the end object, this prevents iterating past the end...
const MeshBase & mesh
Definition: dtk_adapter.h:96
virtual element_iterator local_elements_end()=0
Teuchos::RCP<DataTransferKit::FieldManager<MeshContainerType> > libMesh::DTKAdapter::get_target_coords ( )

Definition at line 66 of file dtk_adapter.h.

References find_sys(), get_element_topology(), get_semi_local_nodes(), get_values_to_fill(), target_coords, and update_variable_values().

Referenced by libMesh::DTKSolutionTransfer::transfer().

66 { return target_coords; }
Teuchos::RCP< DataTransferKit::FieldManager< MeshContainerType > > target_coords
Definition: dtk_adapter.h:104
Teuchos::RCP< DataTransferKit::FieldManager< DTKAdapter::FieldContainerType > > libMesh::DTKAdapter::get_values_to_fill ( std::string  var_name)

Definition at line 192 of file dtk_adapter.C.

References comm, num_local_nodes, and values_to_fill.

Referenced by get_target_coords(), and libMesh::DTKSolutionTransfer::transfer().

193 {
194  if (values_to_fill.find(var_name) == values_to_fill.end())
195  {
196  Teuchos::ArrayRCP<double> data_space(num_local_nodes);
197  Teuchos::RCP<FieldContainerType> field_container = Teuchos::rcp(new FieldContainerType(data_space, 1));
198  values_to_fill[var_name] = Teuchos::rcp(new DataTransferKit::FieldManager<FieldContainerType>(field_container, comm));
199  }
200 
201  return values_to_fill[var_name];
202 }
Teuchos::RCP< const Teuchos::Comm< int > > comm
Definition: dtk_adapter.h:94
unsigned int num_local_nodes
Definition: dtk_adapter.h:99
std::map< std::string, Teuchos::RCP< DataTransferKit::FieldManager< FieldContainerType > > > values_to_fill
Map of variable names to arrays to be filled by a transfer.
Definition: dtk_adapter.h:107
DataTransferKit::FieldContainer< double > FieldContainerType
Definition: dtk_adapter.h:58
DTKAdapter::RCP_Evaluator libMesh::DTKAdapter::get_variable_evaluator ( std::string  var_name)

Definition at line 178 of file dtk_adapter.C.

References evaluators, find_sys(), and libMesh::sys.

Referenced by get_mesh_manager(), and libMesh::DTKSolutionTransfer::transfer().

179 {
180  if (evaluators.find(var_name) == evaluators.end()) // We haven't created an evaluator for the variable yet
181  {
182  System * sys = find_sys(var_name);
183 
184  // Create the FieldEvaluator
185  evaluators[var_name] = Teuchos::rcp(new DTKEvaluator(*sys, var_name));
186  }
187 
188  return evaluators[var_name];
189 }
std::map< std::string, RCP_Evaluator > evaluators
Map of variable names to RCP_Evaluator objects.
Definition: dtk_adapter.h:110
ImplicitSystem & sys
System * find_sys(std::string var_name)
Small helper function for finding the system containing the variable.
Definition: dtk_adapter.C:239
void libMesh::DTKAdapter::update_variable_values ( std::string  var_name)

After computing values for a variable in this EquationSystems we need to take those values and put them in the actual solution vector.

Definition at line 205 of file dtk_adapter.C.

References libMesh::DofObject::dof_number(), find_sys(), mesh, libMesh::MeshBase::node_ref(), libMesh::System::number(), libMesh::ParallelObject::processor_id(), libMesh::DofObject::processor_id(), libMesh::System::solution, libMesh::sys, values_to_fill, libMesh::System::variable_number(), and vertices.

Referenced by get_target_coords(), and libMesh::DTKSolutionTransfer::transfer().

206 {
207  System * sys = find_sys(var_name);
208  unsigned int var_num = sys->variable_number(var_name);
209 
210  Teuchos::RCP<FieldContainerType> values = values_to_fill[var_name]->field();
211 
212  unsigned int i=0;
213  // Loop over the values (one for each node) and assign the value of this variable at each node
214  for (FieldContainerType::iterator it=values->begin(); it != values->end(); ++it)
215  {
216  unsigned int node_num = vertices[i];
217  const Node & node = mesh.node_ref(node_num);
218 
219  if (node.processor_id() == sys->processor_id())
220  {
221  // The 0 is for the component... this only works for LAGRANGE!
222  dof_id_type dof = node.dof_number(sys->number(), var_num, 0);
223  sys->solution->set(dof, *it);
224  }
225 
226  i++;
227  }
228 
229  sys->solution->close();
230 }
ImplicitSystem & sys
virtual const Node & node_ref(const dof_id_type i) const
Definition: mesh_base.h:420
const MeshBase & mesh
Definition: dtk_adapter.h:96
std::map< std::string, Teuchos::RCP< DataTransferKit::FieldManager< FieldContainerType > > > values_to_fill
Map of variable names to arrays to be filled by a transfer.
Definition: dtk_adapter.h:107
Teuchos::ArrayRCP< int > vertices
Definition: dtk_adapter.h:100
System * find_sys(std::string var_name)
Small helper function for finding the system containing the variable.
Definition: dtk_adapter.C:239
uint8_t dof_id_type
Definition: id_types.h:64

Member Data Documentation

Teuchos::RCP<const Teuchos::Comm<int> > libMesh::DTKAdapter::comm
protected

Definition at line 94 of file dtk_adapter.h.

Referenced by DTKAdapter(), and get_values_to_fill().

unsigned int libMesh::DTKAdapter::dim
protected

Definition at line 97 of file dtk_adapter.h.

Referenced by DTKAdapter().

EquationSystems& libMesh::DTKAdapter::es
protected

Definition at line 95 of file dtk_adapter.h.

Referenced by find_sys().

std::map<std::string, RCP_Evaluator> libMesh::DTKAdapter::evaluators
protected

Map of variable names to RCP_Evaluator objects.

Definition at line 110 of file dtk_adapter.h.

Referenced by get_variable_evaluator().

RCP_Evaluator libMesh::DTKAdapter::field_evaluator
protected

Definition at line 103 of file dtk_adapter.h.

const MeshBase& libMesh::DTKAdapter::mesh
protected

Definition at line 96 of file dtk_adapter.h.

Referenced by DTKAdapter(), get_semi_local_nodes(), and update_variable_values().

Teuchos::RCP<DataTransferKit::MeshManager<MeshContainerType> > libMesh::DTKAdapter::mesh_manager
protected

Definition at line 102 of file dtk_adapter.h.

Referenced by DTKAdapter(), and get_mesh_manager().

unsigned int libMesh::DTKAdapter::num_local_nodes
protected

Definition at line 99 of file dtk_adapter.h.

Referenced by DTKAdapter(), and get_values_to_fill().

Teuchos::RCP<DataTransferKit::FieldManager<MeshContainerType> > libMesh::DTKAdapter::target_coords
protected

Definition at line 104 of file dtk_adapter.h.

Referenced by DTKAdapter(), and get_target_coords().

std::map<std::string, Teuchos::RCP<DataTransferKit::FieldManager<FieldContainerType> > > libMesh::DTKAdapter::values_to_fill
protected

Map of variable names to arrays to be filled by a transfer.

Definition at line 107 of file dtk_adapter.h.

Referenced by get_values_to_fill(), and update_variable_values().

Teuchos::ArrayRCP<int> libMesh::DTKAdapter::vertices
protected

Definition at line 100 of file dtk_adapter.h.

Referenced by DTKAdapter(), and update_variable_values().


The documentation for this class was generated from the following files: