libMesh
dtk_adapter.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 #include "libmesh/libmesh_config.h"
21 
22 #ifdef LIBMESH_TRILINOS_HAVE_DTK
23 
24 #include "libmesh/dtk_adapter.h"
25 
26 #include "libmesh/dtk_evaluator.h"
27 #include "libmesh/mesh.h"
28 #include "libmesh/numeric_vector.h"
29 #include "libmesh/elem.h"
30 #include "libmesh/equation_systems.h"
31 
32 #include "libmesh/ignore_warnings.h"
33 #include <DTK_MeshTypes.hpp>
34 #include <Teuchos_Comm.hpp>
35 #include "libmesh/restore_warnings.h"
36 
37 #include <vector>
38 
39 namespace libMesh
40 {
41 
42 DTKAdapter::DTKAdapter(Teuchos::RCP<const Teuchos::Comm<int>> in_comm, EquationSystems & in_es):
43  comm(in_comm),
44  es(in_es),
45  mesh(in_es.get_mesh()),
46  dim(mesh.mesh_dimension())
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 
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 }
176 
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 }
190 
191 Teuchos::RCP<DataTransferKit::FieldManager<DTKAdapter::FieldContainerType>>
192 DTKAdapter::get_values_to_fill(std::string var_name)
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 }
203 
204 void
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 }
231 
232 
238 System *
239 DTKAdapter::find_sys(std::string var_name)
240 {
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 }
257 
258 DataTransferKit::DTK_ElementTopology
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 }
278 
279 void
280 DTKAdapter::get_semi_local_nodes(std::set<unsigned int> & semi_local_nodes)
281 {
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 }
293 
294 } // namespace libMesh
295 
296 #endif // #ifdef LIBMESH_TRILINOS_HAVE_DTK
This is the EquationSystems class.
std::map< std::string, RCP_Evaluator > evaluators
Map of variable names to RCP_Evaluator objects.
Definition: dtk_adapter.h:110
A Node is like a Point, but with more information.
Definition: node.h:52
virtual ElemType type() const =0
Teuchos::RCP< const Teuchos::Comm< int > > comm
Definition: dtk_adapter.h:94
The definition of the const_element_iterator struct.
Definition: mesh_base.h:1494
unsigned int dim
ImplicitSystem & sys
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
ElemType
Defines an enum for geometric element types.
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...
Definition: dtk_adapter.C:205
unsigned int num_local_nodes
Definition: dtk_adapter.h:99
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
MeshBase & mesh
const class libmesh_nullptr_t libmesh_nullptr
dof_id_type dof_number(const unsigned int s, const unsigned int var, const unsigned int comp) const
Definition: dof_object.h:810
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...
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...
Definition: dtk_evaluator.h:60
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
libmesh_assert(j)
virtual unsigned int n_nodes() const =0
EquationSystems & es
Definition: dtk_adapter.h:95
unsigned short int variable_number(const std::string &var) const
Definition: system.C:1263
RCP_Evaluator get_variable_evaluator(std::string var_name)
Definition: dtk_adapter.C:178
This is the base class for classes which contain information related to any physical process that mig...
Definition: system.h:76
virtual element_iterator local_elements_end()=0
UniquePtr< NumericVector< Number > > solution
Data structure to hold solution values.
Definition: system.h:1523
DTKAdapter(Teuchos::RCP< const Teuchos::Comm< int >> in_comm, EquationSystems &in_es)
Definition: dtk_adapter.C:42
Teuchos::RCP< EvaluatorType > RCP_Evaluator
Definition: dtk_adapter.h:61
unsigned int n_systems() const
dof_id_type n_local_elem() const
Definition: mesh_base.h:372
unsigned int number() const
Definition: system.h:2006
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::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
const T_sys & get_system(const std::string &name) const
dof_id_type node_id(const unsigned int i) const
Definition: elem.h:1831
dof_id_type id() const
Definition: dof_object.h:632
Teuchos::RCP< DataTransferKit::MeshManager< MeshContainerType > > mesh_manager
Definition: dtk_adapter.h:102
Teuchos::ArrayRCP< int > vertices
Definition: dtk_adapter.h:100
DataTransferKit::FieldContainer< double > FieldContainerType
Definition: dtk_adapter.h:58
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
processor_id_type processor_id() const
processor_id_type processor_id() const
Definition: dof_object.h:694
Teuchos::RCP< DataTransferKit::FieldManager< FieldContainerType > > get_values_to_fill(std::string var_name)
Definition: dtk_adapter.C:192