exodusII_io_helper.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2018 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 
20 
21 
22 #ifdef LIBMESH_HAVE_EXODUS_API
23 
24 #include <algorithm>
25 #include <functional>
26 #include <sstream>
27 #include <cstdlib> // std::strtol
28 
29 #include "libmesh/boundary_info.h"
30 #include "libmesh/enum_elem_type.h"
31 #include "libmesh/elem.h"
32 #include "libmesh/system.h"
33 #include "libmesh/numeric_vector.h"
34 #include "libmesh/string_to_enum.h"
35 #include "libmesh/enum_elem_type.h"
36 #include "libmesh/int_range.h"
37 
38 #ifdef DEBUG
39 #include "libmesh/mesh_tools.h" // for elem_types warning
40 #endif
41 
42 // This macro returns the length of the array a. Don't
43 // try using it on empty arrays, since it accesses the
44 // zero'th element.
45 #define ARRAY_LENGTH(a) (sizeof((a))/sizeof((a)[0]))
46 
47 // Anonymous namespace for file local data
48 namespace
49 {
50 using namespace libMesh;
51 
52 // Define equivalence classes of Cubit/Exodus element types that map to
53 // libmesh ElemTypes
54 std::map<std::string, ElemType> element_equivalence_map;
55 
56 // This function initializes the element_equivalence_map the first time it
57 // is called, and returns early all other times.
58 void init_element_equivalence_map()
59 {
60  if (element_equivalence_map.empty())
61  {
62  // We use an ExodusII SPHERE element to represent a NodeElem
63  element_equivalence_map["SPHERE"] = NODEELEM;
64 
65  // EDGE2 equivalences
66  element_equivalence_map["EDGE"] = EDGE2;
67  element_equivalence_map["EDGE2"] = EDGE2;
68  element_equivalence_map["TRUSS"] = EDGE2;
69  element_equivalence_map["BEAM"] = EDGE2;
70  element_equivalence_map["BAR"] = EDGE2;
71  element_equivalence_map["TRUSS2"] = EDGE2;
72  element_equivalence_map["BEAM2"] = EDGE2;
73  element_equivalence_map["BAR2"] = EDGE2;
74 
75  // EDGE3 equivalences
76  element_equivalence_map["EDGE3"] = EDGE3;
77  element_equivalence_map["TRUSS3"] = EDGE3;
78  element_equivalence_map["BEAM3"] = EDGE3;
79  element_equivalence_map["BAR3"] = EDGE3;
80 
81  // QUAD4 equivalences
82  element_equivalence_map["QUAD"] = QUAD4;
83  element_equivalence_map["QUAD4"] = QUAD4;
84 
85  // QUADSHELL4 equivalences
86  element_equivalence_map["SHELL"] = QUADSHELL4;
87  element_equivalence_map["SHELL4"] = QUADSHELL4;
88 
89  // QUAD8 equivalences
90  element_equivalence_map["QUAD8"] = QUAD8;
91 
92  // QUADSHELL8 equivalences
93  element_equivalence_map["SHELL8"] = QUADSHELL8;
94 
95  // QUAD9 equivalences
96  element_equivalence_map["QUAD9"] = QUAD9;
97  // element_equivalence_map["SHELL9"] = QUAD9;
98 
99  // TRI3 equivalences
100  element_equivalence_map["TRI"] = TRI3;
101  element_equivalence_map["TRI3"] = TRI3;
102  element_equivalence_map["TRIANGLE"] = TRI3;
103 
104  // TRISHELL3 equivalences
105  element_equivalence_map["TRISHELL"] = TRISHELL3;
106  element_equivalence_map["TRISHELL3"] = TRISHELL3;
107 
108  // TRI6 equivalences
109  element_equivalence_map["TRI6"] = TRI6;
110  // element_equivalence_map["TRISHELL6"] = TRI6;
111 
112  // HEX8 equivalences
113  element_equivalence_map["HEX"] = HEX8;
114  element_equivalence_map["HEX8"] = HEX8;
115 
116  // HEX20 equivalences
117  element_equivalence_map["HEX20"] = HEX20;
118 
119  // HEX27 equivalences
120  element_equivalence_map["HEX27"] = HEX27;
121 
122  // TET4 equivalences
123  element_equivalence_map["TETRA"] = TET4;
124  element_equivalence_map["TETRA4"] = TET4;
125 
126  // TET10 equivalences
127  element_equivalence_map["TETRA10"] = TET10;
128 
129  // PRISM6 equivalences
130  element_equivalence_map["WEDGE"] = PRISM6;
131 
132  // PRISM15 equivalences
133  element_equivalence_map["WEDGE15"] = PRISM15;
134 
135  // PRISM18 equivalences
136  element_equivalence_map["WEDGE18"] = PRISM18;
137 
138  // PYRAMID5 equivalences
139  element_equivalence_map["PYRAMID"] = PYRAMID5;
140  element_equivalence_map["PYRAMID5"] = PYRAMID5;
141 
142  // PYRAMID13 equivalences
143  element_equivalence_map["PYRAMID13"] = PYRAMID13;
144 
145  // PYRAMID14 equivalences
146  element_equivalence_map["PYRAMID14"] = PYRAMID14;
147  }
148 }
149 }
150 
151 
152 
153 namespace libMesh
154 {
155 
156 // ------------------------------------------------------------
157 // ExodusII_IO_Helper::ElementMaps static data
158 
159 // 0D node map definitions
161 
162 // 1D node map definitions
164 const int ExodusII_IO_Helper::ElementMaps::edge3_node_map[3] = {0, 1, 2};
165 
166 // 1D edge maps
167 // FIXME: This notion may or may not be defined in ExodusII
170 
171 // 2D node map definitions
172 const int ExodusII_IO_Helper::ElementMaps::quad4_node_map[4] = {0, 1, 2, 3};
173 const int ExodusII_IO_Helper::ElementMaps::quad8_node_map[8] = {0, 1, 2, 3, 4, 5, 6, 7};
174 const int ExodusII_IO_Helper::ElementMaps::quad9_node_map[9] = {0, 1, 2, 3, 4, 5, 6, 7, 8};
175 const int ExodusII_IO_Helper::ElementMaps::tri3_node_map[3] = {0, 1, 2};
176 const int ExodusII_IO_Helper::ElementMaps::tri6_node_map[6] = {0, 1, 2, 3, 4, 5};
177 
178 // 2D edge map definitions
179 const int ExodusII_IO_Helper::ElementMaps::tri_edge_map[3] = {0, 1, 2};
180 const int ExodusII_IO_Helper::ElementMaps::quad_edge_map[4] = {0, 1, 2, 3};
182 const int ExodusII_IO_Helper::ElementMaps::quadshell4_edge_map[4] = {0, 1, 2, 3};
183 
184 //These take a libMesh ID and turn it into an Exodus ID
189 
190 // 3D node map definitions
191 const int ExodusII_IO_Helper::ElementMaps::hex8_node_map[8] = {0, 1, 2, 3, 4, 5, 6, 7};
192 const int ExodusII_IO_Helper::ElementMaps::hex20_node_map[20] = { 0, 1, 2, 3, 4, 5, 6, 7, 8, 9,
193  10, 11, 12, 13, 14, 15, 16, 17, 18, 19};
194 
195 // Perhaps an older Hex27 node numbering? This no longer works.
196 //const int ExodusII_IO_Helper::ElementMaps::hex27_node_map[27] = { 1, 5, 6, 2, 0, 4, 7, 3, 13, 17, 14, 9, 8, 16,
197 // 18, 10, 12, 19, 15, 11, 24, 25, 22, 26, 21, 23, 20};
198 
199 //DRG: Using the newest exodus documentation available on sourceforge and using Table 2 to see
200 // where all of the nodes over 21 are supposed to go... we come up with:
202  // Vertex and mid-edge nodes
203  0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19,
204  // Mid-face nodes and centroid
205  21, 25, 24, 26, 23, 22, 20};
206 //20 21 22 23 24 25 26 // LibMesh indices
207 
208 // The hex27 appears to be the only element without a 1:1 map between its
209 // node numbering and libmesh's. Therefore when writing out hex27's we need
210 // to invert this map...
212  // Vertex and mid-edge nodes
213  0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19,
214  // Mid-face nodes and centroid
215  26, 20, 25, 24, 22, 21, 23};
216 //20 21 22 23 24 25 26
217 
218 
219 const int ExodusII_IO_Helper::ElementMaps::tet4_node_map[4] = {0, 1, 2, 3};
220 const int ExodusII_IO_Helper::ElementMaps::tet10_node_map[10] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9};
221 
222 const int ExodusII_IO_Helper::ElementMaps::prism6_node_map[6] = {0, 1, 2, 3, 4, 5};
223 const int ExodusII_IO_Helper::ElementMaps::prism15_node_map[15] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9,
224  10, 11, 12, 13, 14};
225 const int ExodusII_IO_Helper::ElementMaps::prism18_node_map[18] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9,
226  10, 11, 12, 13, 14, 15, 16, 17};
227 const int ExodusII_IO_Helper::ElementMaps::pyramid5_node_map[5] = {0, 1, 2, 3, 4};
228 const int ExodusII_IO_Helper::ElementMaps::pyramid13_node_map[13] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12};
229 const int ExodusII_IO_Helper::ElementMaps::pyramid14_node_map[14] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13};
230 
231 // Shell element face maps
234 
235 //These take a libMesh ID and turn it into an Exodus ID
238 
239 // 3D face map definitions
240 const int ExodusII_IO_Helper::ElementMaps::tet_face_map[4] = {1, 2, 3, 0};
241 
242 const int ExodusII_IO_Helper::ElementMaps::hex_face_map[6] = {1, 2, 3, 4, 0, 5};
243 const int ExodusII_IO_Helper::ElementMaps::hex27_face_map[6] = {1, 2, 3, 4, 0, 5};
244 //const int ExodusII_IO_Helper::ElementMaps::hex27_face_map[6] = {1, 0, 3, 5, 4, 2};
245 const int ExodusII_IO_Helper::ElementMaps::prism_face_map[5] = {1, 2, 3, 0, 4};
246 
247 // Pyramids are not included in the ExodusII specification. The ordering below matches
248 // the sideset ordering that CUBIT generates.
249 const int ExodusII_IO_Helper::ElementMaps::pyramid_face_map[5] = {0, 1, 2, 3, 4};
250 
251 //These take a libMesh ID and turn it into an Exodus ID
253 const int ExodusII_IO_Helper::ElementMaps::hex_inverse_face_map[6] = {5, 1, 2, 3, 4, 6};
254 const int ExodusII_IO_Helper::ElementMaps::hex27_inverse_face_map[6] = {5, 1, 2, 3, 4, 6};
255 //const int ExodusII_IO_Helper::ElementMaps::hex27_inverse_face_map[6] = {2, 1, 6, 3, 5, 4};
256 const int ExodusII_IO_Helper::ElementMaps::prism_inverse_face_map[5] = {4, 1, 2, 3, 5};
257 
258 // Pyramids are not included in the ExodusII specification. The ordering below matches
259 // the sideset ordering that CUBIT generates.
260 const int ExodusII_IO_Helper::ElementMaps::pyramid_inverse_face_map[5] = {1, 2, 3, 4, 5};
261 
262 
263 
264 // ExodusII_IO_Helper::Conversion static data
266 
267 // ------------------------------------------------------------
268 // ExodusII_IO_Helper class members
269 
271  bool v,
272  bool run_only_on_proc0,
273  bool single_precision) :
274  ParallelObject(parent),
275  ex_id(0),
276  ex_err(0),
277  num_dim(0),
278  num_global_vars(0),
279  num_nodes(0),
280  num_elem(0),
281  num_elem_blk(0),
282  num_node_sets(0),
283  num_side_sets(0),
284  num_elem_this_blk(0),
285  num_nodes_per_elem(0),
286  num_attr(0),
287  num_elem_all_sidesets(0),
288  num_time_steps(0),
289  num_nodal_vars(0),
290  num_elem_vars(0),
291  verbose(v),
292  opened_for_writing(false),
293  opened_for_reading(false),
294  _run_only_on_proc0(run_only_on_proc0),
295  _elem_vars_initialized(false),
296  _global_vars_initialized(false),
297  _nodal_vars_initialized(false),
298  _use_mesh_dimension_instead_of_spatial_dimension(false),
299  _write_as_dimension(0),
300  _single_precision(single_precision)
301 {
302  title.resize(MAX_LINE_LENGTH+1);
303  elem_type.resize(MAX_STR_LENGTH);
304 }
305 
306 
307 
309 {
310 }
311 
312 
313 
315 {
316  return elem_type.data();
317 }
318 
319 
320 
321 void ExodusII_IO_Helper::message(const std::string & msg)
322 {
323  if (verbose) libMesh::out << msg << std::endl;
324 }
325 
326 
327 
328 void ExodusII_IO_Helper::message(const std::string & msg, int i)
329 {
330  if (verbose) libMesh::out << msg << i << "." << std::endl;
331 }
332 
333 
334 
335 void ExodusII_IO_Helper::open(const char * filename, bool read_only)
336 {
337  // Version of Exodus you are using
338  float ex_version = 0.;
339 
340  int comp_ws = 0;
341 
342  if (_single_precision)
343  comp_ws = cast_int<int>(sizeof(float));
344 
345  // Fall back on double precision when necessary since ExodusII
346  // doesn't seem to support long double
347  else
348  comp_ws = cast_int<int>(std::min(sizeof(Real), sizeof(double)));
349 
350  // Word size in bytes of the floating point data as they are stored
351  // in the ExodusII file. "If this argument is 0, the word size of the
352  // floating point data already stored in the file is returned"
353  int io_ws = 0;
354 
355  ex_id = exII::ex_open(filename,
356  read_only ? EX_READ : EX_WRITE,
357  &comp_ws,
358  &io_ws,
359  &ex_version);
360 
361  std::string err_msg = std::string("Error opening ExodusII mesh file: ") + std::string(filename);
362  EX_CHECK_ERR(ex_id, err_msg);
363  if (verbose) libMesh::out << "File opened successfully." << std::endl;
364 
365  if (read_only)
366  opened_for_reading = true;
367  else
368  opened_for_writing = true;
369 
370  current_filename = std::string(filename);
371 }
372 
373 
374 
376 {
378  title.data(),
379  &num_dim,
380  &num_nodes,
381  &num_elem,
382  &num_elem_blk,
383  &num_node_sets,
384  &num_side_sets);
385 
386  EX_CHECK_ERR(ex_err, "Error retrieving header info.");
387 
388  this->read_num_time_steps();
389 
391  EX_CHECK_ERR(ex_err, "Error reading number of nodal variables.");
392 
394  EX_CHECK_ERR(ex_err, "Error reading number of elemental variables.");
395 
397  EX_CHECK_ERR(ex_err, "Error reading number of global variables.");
398 
399  message("Exodus header info retrieved successfully.");
400 }
401 
402 
403 
404 
406 {
407  // The QA records are four MAX_STR_LENGTH-byte character strings.
408  int num_qa_rec =
409  this->inquire(exII::EX_INQ_QA, "Error retrieving number of QA records");
410 
411  if (verbose)
412  libMesh::out << "Found "
413  << num_qa_rec
414  << " QA record(s) in the Exodus file."
415  << std::endl;
416 
417  if (num_qa_rec > 0)
418  {
419  // How to dynamically allocate an array of fixed-size char * arrays in C++.
420  // http://stackoverflow.com/questions/8529359/creating-a-dynamic-sized-array-of-fixed-sized-int-arrays-in-c
421  typedef char * inner_array_t[4];
422  inner_array_t * qa_record = new inner_array_t[num_qa_rec];
423 
424  for (int i=0; i<num_qa_rec; i++)
425  for (int j=0; j<4; j++)
426  qa_record[i][j] = new char[MAX_STR_LENGTH+1];
427 
428  ex_err = exII::ex_get_qa (ex_id, qa_record);
429  EX_CHECK_ERR(ex_err, "Error reading the QA records.");
430 
431  // Print the QA records
432  if (verbose)
433  {
434  for (int i=0; i<num_qa_rec; i++)
435  {
436  libMesh::out << "QA Record: " << i << std::endl;
437  for (int j=0; j<4; j++)
438  libMesh::out << qa_record[i][j] << std::endl;
439  }
440  }
441 
442 
443  // Clean up dynamically-allocated memory
444  for (int i=0; i<num_qa_rec; i++)
445  for (int j=0; j<4; j++)
446  delete [] qa_record[i][j];
447 
448  delete [] qa_record;
449  }
450 }
451 
452 
453 
454 
456 {
457  if (verbose)
458  libMesh::out << "Title: \t" << title.data() << std::endl
459  << "Mesh Dimension: \t" << num_dim << std::endl
460  << "Number of Nodes: \t" << num_nodes << std::endl
461  << "Number of elements: \t" << num_elem << std::endl
462  << "Number of elt blocks: \t" << num_elem_blk << std::endl
463  << "Number of node sets: \t" << num_node_sets << std::endl
464  << "Number of side sets: \t" << num_side_sets << std::endl;
465 }
466 
467 
468 
470 {
471  x.resize(num_nodes);
472  y.resize(num_nodes);
473  z.resize(num_nodes);
474 
475  if (num_nodes)
476  {
478  static_cast<void *>(x.data()),
479  static_cast<void *>(y.data()),
480  static_cast<void *>(z.data()));
481 
482  EX_CHECK_ERR(ex_err, "Error retrieving nodal data.");
483  message("Nodal data retrieved successfully.");
484  }
485 }
486 
487 
488 
490 {
491  node_num_map.resize(num_nodes);
492 
494  node_num_map.empty() ? nullptr : node_num_map.data());
495 
496  EX_CHECK_ERR(ex_err, "Error retrieving nodal number map.");
497  message("Nodal numbering map retrieved successfully.");
498 
499  if (verbose)
500  {
501  libMesh::out << "[" << this->processor_id() << "] node_num_map[i] = ";
502  for (unsigned int i=0; i<static_cast<unsigned int>(std::min(10, num_nodes-1)); ++i)
503  libMesh::out << node_num_map[i] << ", ";
504  libMesh::out << "... " << node_num_map.back() << std::endl;
505  }
506 }
507 
508 
509 void ExodusII_IO_Helper::print_nodes(std::ostream & out_stream)
510 {
511  for (int i=0; i<num_nodes; i++)
512  out_stream << "(" << x[i] << ", " << y[i] << ", " << z[i] << ")" << std::endl;
513 }
514 
515 
516 
518 {
519  block_ids.resize(num_elem_blk);
520  // Get all element block IDs.
522  block_ids.empty() ? nullptr : block_ids.data());
523  // Usually, there is only one
524  // block since there is only
525  // one type of element.
526  // However, there could be more.
527 
528  EX_CHECK_ERR(ex_err, "Error getting block IDs.");
529  message("All block IDs retrieved successfully.");
530 
531  char name_buffer[MAX_STR_LENGTH+1];
532  for (int i=0; i<num_elem_blk; ++i)
533  {
535  block_ids[i], name_buffer);
536  EX_CHECK_ERR(ex_err, "Error getting block name.");
537  id_to_block_names[block_ids[i]] = name_buffer;
538  }
539  message("All block names retrieved successfully.");
540 }
541 
542 
543 
545 {
546  libmesh_assert_less (static_cast<unsigned int>(index), block_ids.size());
547 
548  return block_ids[index];
549 }
550 
551 
552 
554 {
555  libmesh_assert_less (static_cast<unsigned int>(index), block_ids.size());
556 
557  return id_to_block_names[block_ids[index]];
558 }
559 
560 
561 
563 {
564  libmesh_assert_less (static_cast<unsigned int>(index), ss_ids.size());
565 
566  return ss_ids[index];
567 }
568 
569 
570 
572 {
573  libmesh_assert_less (static_cast<unsigned int>(index), ss_ids.size());
574 
575  return id_to_ss_names[ss_ids[index]];
576 }
577 
578 
579 
581 {
582  libmesh_assert_less (static_cast<unsigned int>(index), nodeset_ids.size());
583 
584  return nodeset_ids[index];
585 }
586 
587 
588 
590 {
591  libmesh_assert_less (static_cast<unsigned int>(index), nodeset_ids.size());
592 
593  return id_to_ns_names[nodeset_ids[index]];
594 }
595 
596 
597 
598 
600 {
601  libmesh_assert_less (static_cast<unsigned int>(block), block_ids.size());
602 
604  block_ids[block],
605  elem_type.data(),
608  &num_attr);
609  if (verbose)
610  libMesh::out << "Reading a block of " << num_elem_this_blk
611  << " " << elem_type.data() << "(s)"
612  << " having " << num_nodes_per_elem
613  << " nodes per element." << std::endl;
614 
615  EX_CHECK_ERR(ex_err, "Error getting block info.");
616  message("Info retrieved successfully for block: ", block);
617 
618 
619 
620  // Read in the connectivity of the elements of this block,
621  // watching out for the case where we actually have no
622  // elements in this block (possible with parallel files)
624 
625  if (!connect.empty())
626  {
628  block_ids[block],
629  connect.data());
630 
631  EX_CHECK_ERR(ex_err, "Error reading block connectivity.");
632  message("Connectivity retrieved successfully for block: ", block);
633  }
634 }
635 
636 
637 
638 
640 {
641  elem_num_map.resize(num_elem);
642 
644  elem_num_map.empty() ? nullptr : elem_num_map.data());
645 
646  EX_CHECK_ERR(ex_err, "Error retrieving element number map.");
647  message("Element numbering map retrieved successfully.");
648 
649 
650  if (verbose)
651  {
652  libMesh::out << "[" << this->processor_id() << "] elem_num_map[i] = ";
653  for (unsigned int i=0; i<static_cast<unsigned int>(std::min(10, num_elem-1)); ++i)
654  libMesh::out << elem_num_map[i] << ", ";
655  libMesh::out << "... " << elem_num_map.back() << std::endl;
656  }
657 }
658 
659 
660 
662 {
663  ss_ids.resize(num_side_sets);
664  if (num_side_sets > 0)
665  {
667  ss_ids.data());
668  EX_CHECK_ERR(ex_err, "Error retrieving sideset information.");
669  message("All sideset information retrieved successfully.");
670 
671  // Resize appropriate data structures -- only do this once outside the loop
674 
675  // Inquire about the length of the concatenated side sets element list
676  num_elem_all_sidesets = inquire(exII::EX_INQ_SS_ELEM_LEN, "Error retrieving length of the concatenated side sets element list!");
677 
681  }
682 
683  char name_buffer[MAX_STR_LENGTH+1];
684  for (int i=0; i<num_side_sets; ++i)
685  {
687  ss_ids[i], name_buffer);
688  EX_CHECK_ERR(ex_err, "Error getting side set name.");
689  id_to_ss_names[ss_ids[i]] = name_buffer;
690  }
691  message("All side set names retrieved successfully.");
692 }
693 
695 {
696  nodeset_ids.resize(num_node_sets);
697  if (num_node_sets > 0)
698  {
700  nodeset_ids.data());
701  EX_CHECK_ERR(ex_err, "Error retrieving nodeset information.");
702  message("All nodeset information retrieved successfully.");
703 
704  // Resize appropriate data structures -- only do this once outnode the loop
707  }
708 
709  char name_buffer[MAX_STR_LENGTH+1];
710  for (int i=0; i<num_node_sets; ++i)
711  {
713  nodeset_ids[i], name_buffer);
714  EX_CHECK_ERR(ex_err, "Error getting node set name.");
715  id_to_ns_names[nodeset_ids[i]] = name_buffer;
716  }
717  message("All node set names retrieved successfully.");
718 }
719 
720 
721 
722 void ExodusII_IO_Helper::read_sideset(int id, int offset)
723 {
724  libmesh_assert_less (static_cast<unsigned int>(id), ss_ids.size());
725  libmesh_assert_less (static_cast<unsigned int>(id), num_sides_per_set.size());
726  libmesh_assert_less (static_cast<unsigned int>(id), num_df_per_set.size());
727  libmesh_assert_less_equal (static_cast<unsigned int>(offset), elem_list.size());
728  libmesh_assert_less_equal (static_cast<unsigned int>(offset), side_list.size());
729 
731  ss_ids[id],
732  &num_sides_per_set[id],
733  &num_df_per_set[id]);
734  EX_CHECK_ERR(ex_err, "Error retrieving sideset parameters.");
735  message("Parameters retrieved successfully for sideset: ", id);
736 
737 
738  // It's OK for offset==elem_list.size() as long as num_sides_per_set[id]==0
739  // because in that case we don't actually read anything...
740 #ifdef DEBUG
741  if (static_cast<unsigned int>(offset) == elem_list.size() ||
742  static_cast<unsigned int>(offset) == side_list.size() )
743  libmesh_assert_equal_to (num_sides_per_set[id], 0);
744 #endif
745 
746 
747  // Don't call ex_get_side_set unless there are actually sides there to get.
748  // Exodus prints an annoying warning in DEBUG mode otherwise...
749  if (num_sides_per_set[id] > 0)
750  {
752  ss_ids[id],
753  &elem_list[offset],
754  &side_list[offset]);
755  EX_CHECK_ERR(ex_err, "Error retrieving sideset data.");
756  message("Data retrieved successfully for sideset: ", id);
757 
758  for (int i=0; i<num_sides_per_set[id]; i++)
759  id_list[i+offset] = ss_ids[id];
760  }
761 }
762 
763 
764 
766 {
767  libmesh_assert_less (static_cast<unsigned int>(id), nodeset_ids.size());
768  libmesh_assert_less (static_cast<unsigned int>(id), num_nodes_per_set.size());
769  libmesh_assert_less (static_cast<unsigned int>(id), num_node_df_per_set.size());
770 
772  nodeset_ids[id],
773  &num_nodes_per_set[id],
774  &num_node_df_per_set[id]);
775  EX_CHECK_ERR(ex_err, "Error retrieving nodeset parameters.");
776  message("Parameters retrieved successfully for nodeset: ", id);
777 
778  node_list.resize(num_nodes_per_set[id]);
779 
780  // Don't call ex_get_node_set unless there are actually nodes there to get.
781  // Exodus prints an annoying warning message in DEBUG mode otherwise...
782  if (num_nodes_per_set[id] > 0)
783  {
785  nodeset_ids[id],
786  node_list.data());
787 
788  EX_CHECK_ERR(ex_err, "Error retrieving nodeset data.");
789  message("Data retrieved successfully for nodeset: ", id);
790  }
791 }
792 
793 
794 
796 {
797  // Always call close on processor 0.
798  // If we're running on multiple processors, i.e. as one of several Nemesis files,
799  // we call close on all processors...
800  if ((this->processor_id() == 0) || (!_run_only_on_proc0))
801  {
802  // Don't close the file if it was never opened, this raises an Exodus error
804  {
806  EX_CHECK_ERR(ex_err, "Error closing Exodus file.");
807  message("Exodus file closed successfully.");
808  }
809  }
810 }
811 
812 
813 
814 int ExodusII_IO_Helper::inquire(int req_info_in, std::string error_msg)
815 {
816  int ret_int = 0;
817  char ret_char = 0;
818  float ret_float = 0.;
819 
821  req_info_in,
822  &ret_int,
823  &ret_float,
824  &ret_char);
825 
826  EX_CHECK_ERR(ex_err, error_msg);
827 
828  return ret_int;
829 }
830 
831 
832 
834 {
835  // Make sure we have an up-to-date count of the number of time steps in the file.
836  this->read_num_time_steps();
837 
838  if (num_time_steps > 0)
839  {
840  time_steps.resize(num_time_steps);
842  EX_CHECK_ERR(ex_err, "Error reading timesteps!");
843  }
844 }
845 
846 
847 
849 {
851  this->inquire(exII::EX_INQ_TIME, "Error retrieving number of time steps");
852 }
853 
854 
855 
856 void ExodusII_IO_Helper::read_nodal_var_values(std::string nodal_var_name, int time_step)
857 {
858  // Read the nodal variable names from file, so we can see if we have the one we're looking for
859  this->read_var_names(NODAL);
860 
861  // See if we can find the variable we are looking for
862  unsigned int var_index = 0;
863  bool found = false;
864 
865  // Do a linear search for nodal_var_name in nodal_var_names
866  for (; var_index<nodal_var_names.size(); ++var_index)
867  {
868  found = (nodal_var_names[var_index] == nodal_var_name);
869  if (found)
870  break;
871  }
872 
873  if (!found)
874  {
875  libMesh::err << "Available variables: " << std::endl;
876  for (const auto & var_name : nodal_var_names)
877  libMesh::err << var_name << std::endl;
878 
879  libmesh_error_msg("Unable to locate variable named: " << nodal_var_name);
880  }
881 
882  // Allocate enough space to store the nodal variable values
883  nodal_var_values.resize(num_nodes);
884 
885  // Call the Exodus API to read the nodal variable values
887  time_step,
888  var_index+1,
889  num_nodes,
890  nodal_var_values.data());
891  EX_CHECK_ERR(ex_err, "Error reading nodal variable values!");
892 }
893 
894 
895 
897 {
898  switch (type)
899  {
900  case NODAL:
902  break;
903  case ELEMENTAL:
905  break;
906  case GLOBAL:
908  break;
909  default:
910  libmesh_error_msg("Unrecognized ExodusVarType " << type);
911  }
912 }
913 
914 
915 
916 void ExodusII_IO_Helper::read_var_names_impl(const char * var_type,
917  int & count,
918  std::vector<std::string> & result)
919 {
920  // First read and store the number of names we have
921  ex_err = exII::ex_get_var_param(ex_id, var_type, &count);
922  EX_CHECK_ERR(ex_err, "Error reading number of variables.");
923 
924  // Do nothing if no variables are detected
925  if (count == 0)
926  return;
927 
928  // Second read the actual names and convert them into a format we can use
929  NamesData names_table(count, MAX_STR_LENGTH);
930 
932  var_type,
933  count,
934  names_table.get_char_star_star()
935  );
936  EX_CHECK_ERR(ex_err, "Error reading variable names!");
937 
938  if (verbose)
939  {
940  libMesh::out << "Read the variable(s) from the file:" << std::endl;
941  for (int i=0; i<count; i++)
942  libMesh::out << names_table.get_char_star(i) << std::endl;
943  }
944 
945  // Allocate enough space for our variable name strings.
946  result.resize(count);
947 
948  // Copy the char buffers into strings.
949  for (int i=0; i<count; i++)
950  result[i] = names_table.get_char_star(i); // calls string::op=(const char *)
951 }
952 
953 
954 
955 
956 void ExodusII_IO_Helper::write_var_names(ExodusVarType type, std::vector<std::string> & names)
957 {
958  switch (type)
959  {
960  case NODAL:
961  this->write_var_names_impl("n", num_nodal_vars, names);
962  break;
963  case ELEMENTAL:
964  this->write_var_names_impl("e", num_elem_vars, names);
965  break;
966  case GLOBAL:
967  this->write_var_names_impl("g", num_global_vars, names);
968  break;
969  default:
970  libmesh_error_msg("Unrecognized ExodusVarType " << type);
971  }
972 }
973 
974 
975 
976 void ExodusII_IO_Helper::write_var_names_impl(const char * var_type, int & count, std::vector<std::string> & names)
977 {
978  // Update the count variable so that it's available to other parts of the class.
979  count = cast_int<int>(names.size());
980 
981  // Write that number of variables to the file.
982  ex_err = exII::ex_put_var_param(ex_id, var_type, count);
983  EX_CHECK_ERR(ex_err, "Error setting number of vars.");
984 
985  if (count > 0)
986  {
987  NamesData names_table(count, MAX_STR_LENGTH);
988 
989  // Store the input names in the format required by Exodus.
990  for (int i=0; i != count; ++i)
991  names_table.push_back_entry(names[i]);
992 
993  if (verbose)
994  {
995  libMesh::out << "Writing variable name(s) to file: " << std::endl;
996  for (int i=0; i != count; ++i)
997  libMesh::out << names_table.get_char_star(i) << std::endl;
998  }
999 
1001  var_type,
1002  count,
1003  names_table.get_char_star_star()
1004  );
1005 
1006  EX_CHECK_ERR(ex_err, "Error writing variable names.");
1007  }
1008 }
1009 
1010 
1011 
1012 void ExodusII_IO_Helper::read_elemental_var_values(std::string elemental_var_name,
1013  int time_step,
1014  std::map<dof_id_type, Real> & elem_var_value_map)
1015 {
1016  this->read_var_names(ELEMENTAL);
1017 
1018  // See if we can find the variable we are looking for
1019  unsigned int var_index = 0;
1020  bool found = false;
1021 
1022  // Do a linear search for nodal_var_name in nodal_var_names
1023  for (; var_index != elem_var_names.size(); ++var_index)
1024  if (elem_var_names[var_index] == elemental_var_name)
1025  {
1026  found = true;
1027  break;
1028  }
1029 
1030  if (!found)
1031  {
1032  libMesh::err << "Available variables: " << std::endl;
1033  for (const auto & var_name : elem_var_names)
1034  libMesh::err << var_name << std::endl;
1035 
1036  libmesh_error_msg("Unable to locate variable named: " << elemental_var_name);
1037  }
1038 
1039  // Sequential index which we can use to look up the element ID in the elem_num_map.
1040  unsigned ex_el_num = 0;
1041 
1042  for (unsigned i=0; i<static_cast<unsigned>(num_elem_blk); i++)
1043  {
1045  block_ids[i],
1046  nullptr,
1048  nullptr,
1049  nullptr);
1050  EX_CHECK_ERR(ex_err, "Error getting number of elements in block.");
1051 
1052  std::vector<Real> block_elem_var_values(num_elem);
1054  time_step,
1055  var_index+1,
1056  block_ids[i],
1058  block_elem_var_values.data());
1059  EX_CHECK_ERR(ex_err, "Error getting elemental values.");
1060 
1061  for (unsigned j=0; j<static_cast<unsigned>(num_elem_this_blk); j++)
1062  {
1063  // Use the elem_num_map to obtain the ID of this element in the Exodus file,
1064  // and remember to subtract 1 since libmesh is zero-based and Exodus is 1-based.
1065  unsigned mapped_elem_id = this->elem_num_map[ex_el_num] - 1;
1066 
1067  // Store the elemental value in the map.
1068  elem_var_value_map[mapped_elem_id] = block_elem_var_values[j];
1069 
1070  // Go to the next sequential element ID.
1071  ex_el_num++;
1072  }
1073  }
1074 }
1075 
1076 
1077 // For Writing Solutions
1078 
1079 void ExodusII_IO_Helper::create(std::string filename)
1080 {
1081  // If we're processor 0, always create the file.
1082  // If we running on all procs, e.g. as one of several Nemesis files, also
1083  // call create there.
1084  if ((this->processor_id() == 0) || (!_run_only_on_proc0))
1085  {
1086  int
1087  comp_ws = 0,
1088  io_ws = 0;
1089 
1090  if (_single_precision)
1091  {
1092  comp_ws = cast_int<int>(sizeof(float));
1093  io_ws = cast_int<int>(sizeof(float));
1094  }
1095  // Fall back on double precision when necessary since ExodusII
1096  // doesn't seem to support long double
1097  else
1098  {
1099  comp_ws = cast_int<int>
1100  (std::min(sizeof(Real), sizeof(double)));
1101  io_ws = cast_int<int>
1102  (std::min(sizeof(Real), sizeof(double)));
1103  }
1104 
1105  ex_id = exII::ex_create(filename.c_str(), EX_CLOBBER, &comp_ws, &io_ws);
1106 
1107  EX_CHECK_ERR(ex_id, "Error creating ExodusII mesh file.");
1108 
1109  if (verbose)
1110  libMesh::out << "File created successfully." << std::endl;
1111  }
1112 
1113  opened_for_writing = true;
1114  current_filename = filename;
1115 }
1116 
1117 
1118 
1119 void ExodusII_IO_Helper::initialize(std::string str_title, const MeshBase & mesh, bool use_discontinuous)
1120 {
1121  // n_active_elem() is a parallel_only function
1122  unsigned int n_active_elem = mesh.n_active_elem();
1123 
1124  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1125  return;
1126 
1127  // If _write_as_dimension is nonzero, use it to set num_dim in the Exodus file.
1128  if (_write_as_dimension)
1131  num_dim = mesh.mesh_dimension();
1132  else
1134 
1135  num_elem = mesh.n_elem();
1136 
1137  if (!use_discontinuous)
1138  {
1139  // Don't rely on mesh.n_nodes() here. If nodes have been
1140  // deleted recently, it will be incorrect.
1141  num_nodes = cast_int<int>(std::distance(mesh.nodes_begin(),
1142  mesh.nodes_end()));
1143  }
1144  else
1145  {
1146  for (const auto & elem : mesh.active_element_ptr_range())
1147  num_nodes += elem->n_nodes();
1148  }
1149 
1150  std::vector<boundary_id_type> unique_side_boundaries;
1151  std::vector<boundary_id_type> unique_node_boundaries;
1152 
1153  mesh.get_boundary_info().build_side_boundary_ids(unique_side_boundaries);
1154  {
1155  // Add shell face boundaries to the list of side boundaries, since ExodusII
1156  // treats these the same way.
1157  std::vector<boundary_id_type> shellface_boundaries;
1158  mesh.get_boundary_info().build_shellface_boundary_ids(shellface_boundaries);
1159  for (const auto & id : shellface_boundaries)
1160  unique_side_boundaries.push_back(id);
1161  }
1162  mesh.get_boundary_info().build_node_boundary_ids(unique_node_boundaries);
1163 
1164  num_side_sets = cast_int<int>(unique_side_boundaries.size());
1165  num_node_sets = cast_int<int>(unique_node_boundaries.size());
1166 
1167  //loop through element and map between block and element vector
1168  std::map<subdomain_id_type, std::vector<unsigned int>> subdomain_map;
1169 
1170  for (const auto & elem : mesh.active_element_ptr_range())
1171  {
1172  // We skip writing infinite elements to the Exodus file, so
1173  // don't put them in the subdomain_map. That way the number of
1174  // blocks should be correct.
1175 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1176  if (elem->infinite())
1177  continue;
1178 #endif
1179 
1180  subdomain_id_type cur_subdomain = elem->subdomain_id();
1181  subdomain_map[cur_subdomain].push_back(elem->id());
1182  }
1183  num_elem_blk = cast_int<int>(subdomain_map.size());
1184 
1185  if (str_title.size() > MAX_LINE_LENGTH)
1186  {
1187  libMesh::err << "Warning, Exodus files cannot have titles longer than "
1188  << MAX_LINE_LENGTH
1189  << " characters. Your title will be truncated."
1190  << std::endl;
1191  str_title.resize(MAX_LINE_LENGTH);
1192  }
1193 
1195  str_title.c_str(),
1196  num_dim,
1197  num_nodes,
1198  n_active_elem,
1199  num_elem_blk,
1200  num_node_sets,
1201  num_side_sets);
1202 
1203  EX_CHECK_ERR(ex_err, "Error initializing new Exodus file.");
1204 }
1205 
1206 
1207 
1208 void ExodusII_IO_Helper::write_nodal_coordinates(const MeshBase & mesh, bool use_discontinuous)
1209 {
1210  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1211  return;
1212 
1213  // Clear existing data from any previous calls.
1214  x.clear();
1215  y.clear();
1216  z.clear();
1217  node_num_map.clear();
1218 
1219  // Reserve space in the nodal coordinate vectors. num_nodes is
1220  // exact, this just allows us to do away with one potentially
1221  // error-inducing loop index.
1222  x.reserve(num_nodes);
1223  y.reserve(num_nodes);
1224  z.reserve(num_nodes);
1225 
1226  // And in the node_num_map - since the nodes aren't organized in
1227  // blocks, libmesh will always write out the identity map
1228  // here... unless there has been some refinement and coarsening, or
1229  // node deletion without a corresponding call to contract(). You
1230  // need to write this any time there could be 'holes' in the node
1231  // numbering, so we write it every time.
1232  node_num_map.reserve(num_nodes);
1233 
1234  // Clear out any previously-mapped node IDs.
1236 
1237  if (!use_discontinuous)
1238  {
1239  for (const auto & node_ptr : mesh.node_ptr_range())
1240  {
1241  const Node & node = *node_ptr;
1242 
1243  x.push_back(node(0) + _coordinate_offset(0));
1244 
1245 #if LIBMESH_DIM > 1
1246  y.push_back(node(1) + _coordinate_offset(1));
1247 #else
1248  y.push_back(0.);
1249 #endif
1250 #if LIBMESH_DIM > 2
1251  z.push_back(node(2) + _coordinate_offset(2));
1252 #else
1253  z.push_back(0.);
1254 #endif
1255 
1256  // Fill in node_num_map entry with the proper (1-based) node id
1257  node_num_map.push_back(node.id() + 1);
1258 
1259  // Also map the zero-based libmesh node id to the 1-based
1260  // Exodus ID it will be assigned (this is equivalent to the
1261  // current size of the x vector).
1262  libmesh_node_num_to_exodus[ cast_int<int>(node.id()) ] = cast_int<int>(x.size());
1263  }
1264  }
1265  else
1266  {
1267  for (const auto & elem : mesh.active_element_ptr_range())
1268  for (unsigned int n=0; n<elem->n_nodes(); n++)
1269  {
1270  x.push_back(elem->point(n)(0));
1271 #if LIBMESH_DIM > 1
1272  y.push_back(elem->point(n)(1));
1273 #else
1274  y.push_back(0.);
1275 #endif
1276 #if LIBMESH_DIM > 2
1277  z.push_back(elem->point(n)(2));
1278 #else
1279  z.push_back(0.);
1280 #endif
1281 
1282  // Let's skip the node_num_map in the discontinuous
1283  // case, since we're effectively duplicating nodes for
1284  // the sake of discontinuous visualization, so it isn't
1285  // clear how to deal with node_num_map here. This means
1286  // that writing discontinuous meshes won't work with
1287  // element numberings that have "holes".
1288  }
1289  }
1290 
1291  if (_single_precision)
1292  {
1293  std::vector<float>
1294  x_single(x.begin(), x.end()),
1295  y_single(y.begin(), y.end()),
1296  z_single(z.begin(), z.end());
1297 
1299  x_single.empty() ? nullptr : x_single.data(),
1300  y_single.empty() ? nullptr : y_single.data(),
1301  z_single.empty() ? nullptr : z_single.data());
1302  }
1303  else
1304  {
1306  x.empty() ? nullptr : x.data(),
1307  y.empty() ? nullptr : y.data(),
1308  z.empty() ? nullptr : z.data());
1309  }
1310 
1311 
1312  EX_CHECK_ERR(ex_err, "Error writing coordinates to Exodus file.");
1313 
1314  if (!use_discontinuous)
1315  {
1316  // Also write the (1-based) node_num_map to the file.
1318  EX_CHECK_ERR(ex_err, "Error writing node_num_map");
1319  }
1320 }
1321 
1322 
1323 
1324 void ExodusII_IO_Helper::write_elements(const MeshBase & mesh, bool use_discontinuous)
1325 {
1326  // n_active_elem() is a parallel_only function
1327  unsigned int n_active_elem = mesh.n_active_elem();
1328 
1329  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1330  return;
1331 
1332  // Map from block ID to a vector of element IDs in that block. Element
1333  // IDs are now of type dof_id_type, subdomain IDs are of type subdomain_id_type.
1334  typedef std::map<subdomain_id_type, std::vector<dof_id_type>> subdomain_map_type;
1335  subdomain_map_type subdomain_map;
1336 
1337  // Loop through element and map between block and element vector.
1338  for (const auto & elem : mesh.active_element_ptr_range())
1339  {
1340  // We skip writing infinite elements to the Exodus file, so
1341  // don't put them in the subdomain_map. That way the number of
1342  // blocks should be correct.
1343 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1344  if (elem->infinite())
1345  continue;
1346 #endif
1347 
1348  subdomain_map[ elem->subdomain_id() ].push_back(elem->id());
1349  }
1350 
1351  // element map vector
1352  num_elem_blk = cast_int<int>(subdomain_map.size());
1353  block_ids.resize(num_elem_blk);
1354  elem_num_map.resize(n_active_elem);
1355  std::vector<int>::iterator curr_elem_map_end = elem_num_map.begin();
1356 
1357  std::vector<int> elem_blk_id;
1358  std::vector<int> num_elem_this_blk_vec;
1359  std::vector<int> num_nodes_per_elem_vec;
1360  std::vector<int> num_attr_vec;
1361  NamesData elem_type_table(num_elem_blk, MAX_STR_LENGTH);
1362 
1363  // Note: It appears that there is a bug in exodusII::ex_put_name where
1364  // the index returned from the ex_id_lkup is erroneously used. For now
1365  // the work around is to use the alternative function ex_put_names, but
1366  // this function requires a char ** data structure.
1367  NamesData names_table(num_elem_blk, MAX_STR_LENGTH);
1368 
1369  // counter indexes into the block_ids vector
1370  unsigned int counter = 0;
1371  for (auto & pr : subdomain_map)
1372  {
1373  block_ids[counter] = pr.first;
1374  names_table.push_back_entry(mesh.subdomain_name(pr.first));
1375 
1376  // Get a reference to a vector of element IDs for this subdomain.
1377  subdomain_map_type::mapped_type & tmp_vec = pr.second;
1378 
1379  // Use the first element in this block to get representative information.
1380  // Note that Exodus assumes all elements in a block are of the same type!
1381  // We are using that same assumption here!
1383  const ExodusII_IO_Helper::Conversion conv =
1384  em.assign_conversion(mesh.elem_ref(tmp_vec[0]).type());
1385  num_nodes_per_elem = mesh.elem_ref(tmp_vec[0]).n_nodes();
1386 
1387  elem_blk_id.push_back(pr.first);
1388  elem_type_table.push_back_entry(conv.exodus_elem_type().c_str());
1389  num_elem_this_blk_vec.push_back(cast_int<int>(tmp_vec.size()));
1390  num_nodes_per_elem_vec.push_back(num_nodes_per_elem);
1391  num_attr_vec.push_back(0); // we don't currently use elem block attributes.
1392  ++counter;
1393  }
1394 
1395  // The "define_maps" parameter should be 0 if node_number_map and
1396  // elem_number_map will not be written later, and nonzero otherwise.
1398  elem_blk_id.data(),
1399  elem_type_table.get_char_star_star(),
1400  num_elem_this_blk_vec.data(),
1401  num_nodes_per_elem_vec.data(),
1402  num_attr_vec.data(),
1403  /*define_maps=*/0);
1404  EX_CHECK_ERR(ex_err, "Error writing element blocks.");
1405 
1406  // This counter is used to fill up the libmesh_elem_num_to_exodus map in the loop below.
1407  unsigned libmesh_elem_num_to_exodus_counter = 0;
1408 
1409  // In the case of discontinuous plotting we initialize a map from (element,node) pairs
1410  // to the corresponding discontinuous node index. This ordering must match the ordering
1411  // used in write_nodal_coordinates.
1412  std::map< std::pair<dof_id_type,unsigned int>, dof_id_type > discontinuous_node_indices;
1413  if (use_discontinuous)
1414  {
1415  dof_id_type node_counter = 1; // Exodus numbering is 1-based
1416  for (const auto & elem : mesh.active_element_ptr_range())
1417  for (unsigned int n=0; n<elem->n_nodes(); n++)
1418  {
1419  std::pair<dof_id_type,unsigned int> id_pair;
1420  id_pair.first = elem->id();
1421  id_pair.second = n;
1422  discontinuous_node_indices[id_pair] = node_counter;
1423  node_counter++;
1424  }
1425  }
1426 
1427  for (auto & pr : subdomain_map)
1428  {
1429  // Get a reference to a vector of element IDs for this subdomain.
1430  subdomain_map_type::mapped_type & tmp_vec = pr.second;
1431 
1432  //Use the first element in this block to get representative information.
1433  //Note that Exodus assumes all elements in a block are of the same type!
1434  //We are using that same assumption here!
1436  const ExodusII_IO_Helper::Conversion conv =
1437  em.assign_conversion(mesh.elem_ref(tmp_vec[0]).type());
1438  num_nodes_per_elem = mesh.elem_ref(tmp_vec[0]).n_nodes();
1439 
1440  connect.resize(tmp_vec.size()*num_nodes_per_elem);
1441 
1442  for (auto i : index_range(tmp_vec))
1443  {
1444  unsigned int elem_id = tmp_vec[i];
1445  libmesh_elem_num_to_exodus[elem_id] = ++libmesh_elem_num_to_exodus_counter; // 1-based indexing for Exodus
1446 
1447  const Elem & elem = mesh.elem_ref(elem_id);
1448 
1449  // We *might* be able to get away with writing mixed element
1450  // types which happen to have the same number of nodes, but
1451  // do we actually *want* to get away with that?
1452  // .) No visualization software would be able to handle it.
1453  // .) There'd be no way for us to read it back in reliably.
1454  // .) Even elements with the same number of nodes may have different connectivities (?)
1455 
1456  // This needs to be more than an assert so we don't fail
1457  // with a mysterious segfault while trying to write mixed
1458  // element meshes in optimized mode.
1459  if (elem.type() != conv.get_canonical_type())
1460  libmesh_error_msg("Error: Exodus requires all elements with a given subdomain ID to be the same type.\n" \
1461  << "Can't write both " \
1462  << Utility::enum_to_string(elem.type()) \
1463  << " and " \
1465  << " in the same block!");
1466 
1467 
1468  for (unsigned int j=0; j<static_cast<unsigned int>(num_nodes_per_elem); ++j)
1469  {
1470  unsigned int connect_index = cast_int<unsigned int>((i*num_nodes_per_elem)+j);
1471  unsigned elem_node_index = conv.get_inverse_node_map(j); // inverse node map is for writing.
1472  if (verbose)
1473  {
1474  libMesh::out << "Exodus node index " << j
1475  << " = LibMesh node index " << elem_node_index << std::endl;
1476  }
1477 
1478  if (!use_discontinuous)
1479  {
1480  // The global id for the current node in libmesh.
1481  dof_id_type libmesh_node_id = elem.node_id(elem_node_index);
1482 
1483  // Find the zero-based libmesh id in the map, this
1484  // should be faster than doing linear searches on
1485  // the node_num_map.
1486  std::map<int, int>::iterator pos =
1487  libmesh_node_num_to_exodus.find(cast_int<int>(libmesh_node_id));
1488 
1489  // Make sure it was found.
1490  if (pos == libmesh_node_num_to_exodus.end())
1491  libmesh_error_msg("libmesh node id " << libmesh_node_id << " not found in node_num_map.");
1492 
1493  // Write the Exodus global node id associated with
1494  // this libmesh node number to the connectivity
1495  // array.
1496  connect[connect_index] = pos->second;
1497  }
1498  else
1499  {
1500  std::pair<dof_id_type,unsigned int> id_pair;
1501  id_pair.first = elem_id;
1502  id_pair.second = elem_node_index;
1503  auto node_it = discontinuous_node_indices.find(id_pair);
1504 
1505  libmesh_assert(node_it != discontinuous_node_indices.end());
1506 
1507  connect[connect_index] = node_it->second;
1508  }
1509  }
1510  }
1511 
1512  ex_err = exII::ex_put_elem_conn(ex_id, pr.first, connect.data());
1513  EX_CHECK_ERR(ex_err, "Error writing element connectivities");
1514 
1515  // This transform command stores its result in a range that begins at the third argument,
1516  // so this command is adding values to the elem_num_map vector starting from curr_elem_map_end.
1517  curr_elem_map_end = std::transform(tmp_vec.begin(),
1518  tmp_vec.end(),
1519  curr_elem_map_end,
1520  std::bind2nd(std::plus<subdomain_map_type::mapped_type::value_type>(), 1)); // Adds one to each id to make a 1-based exodus file!
1521 
1522  // But if we don't want to add one, we just want to put the values
1523  // of tmp_vec into elem_map in the right location, we can use
1524  // std::copy().
1525  // curr_elem_map_end = std::copy(tmp_vec.begin(), tmp_vec.end(), curr_elem_map_end);
1526  }
1527 
1528  // write out the element number map that we created
1530  EX_CHECK_ERR(ex_err, "Error writing element map");
1531 
1532  // Write out the block names
1533  if (num_elem_blk > 0)
1534  {
1536  EX_CHECK_ERR(ex_err, "Error writing element names");
1537  }
1538 
1539 }
1540 
1541 
1542 
1543 
1545 {
1546  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1547  return;
1548 
1550 
1551  // Maps from sideset id to the element and sides
1552  std::map<int, std::vector<int>> elem;
1553  std::map<int, std::vector<int>> side;
1554  std::vector<boundary_id_type> side_boundary_ids;
1555 
1556  {
1557  // Accumulate the vectors to pass into ex_put_side_set
1558  // build_side_list() returns a vector of (elem, side, bc) tuples.
1559  for (const auto & t : mesh.get_boundary_info().build_side_list())
1560  {
1561  std::vector<const Elem *> family;
1562 #ifdef LIBMESH_ENABLE_AMR
1563 
1567  mesh.elem_ref(std::get<0>(t)).active_family_tree_by_side(family, std::get<1>(t), false);
1568 #else
1569  family.push_back(mesh.elem_ptr(std::get<0>(t)));
1570 #endif
1571 
1572  for (const auto & f : family)
1573  {
1574  const ExodusII_IO_Helper::Conversion conv =
1575  em.assign_conversion(mesh.elem_ptr(f->id())->type());
1576 
1577  // Use the libmesh to exodus data structure map to get the proper sideset IDs
1578  // The data structure contains the "collapsed" contiguous ids
1579  elem[std::get<2>(t)].push_back(libmesh_elem_num_to_exodus[f->id()]);
1580  side[std::get<2>(t)].push_back(conv.get_inverse_side_map(std::get<1>(t)));
1581  }
1582  }
1583 
1584  mesh.get_boundary_info().build_side_boundary_ids(side_boundary_ids);
1585  }
1586 
1587  {
1588  // add data for shell faces, if needed
1589 
1590  // Accumulate the vectors to pass into ex_put_side_set
1591  for (const auto & t : mesh.get_boundary_info().build_shellface_list())
1592  {
1593  std::vector<const Elem *> family;
1594 #ifdef LIBMESH_ENABLE_AMR
1595 
1599  mesh.elem_ref(std::get<0>(t)).active_family_tree_by_side(family, std::get<1>(t), false);
1600 #else
1601  family.push_back(mesh.elem_ptr(std::get<0>(t)));
1602 #endif
1603 
1604  for (const auto & f : family)
1605  {
1606  const ExodusII_IO_Helper::Conversion conv =
1607  em.assign_conversion(mesh.elem_ptr(f->id())->type());
1608 
1609  // Use the libmesh to exodus data structure map to get the proper sideset IDs
1610  // The data structure contains the "collapsed" contiguous ids
1611  elem[std::get<2>(t)].push_back(libmesh_elem_num_to_exodus[f->id()]);
1612  side[std::get<2>(t)].push_back(conv.get_inverse_shellface_map(std::get<1>(t)));
1613  }
1614  }
1615 
1616  std::vector<boundary_id_type> shellface_boundary_ids;
1617  mesh.get_boundary_info().build_shellface_boundary_ids(shellface_boundary_ids);
1618  for (const auto & id : shellface_boundary_ids)
1619  side_boundary_ids.push_back(id);
1620  }
1621 
1622  // Write out the sideset names, but only if there is something to write
1623  if (side_boundary_ids.size() > 0)
1624  {
1625  NamesData names_table(side_boundary_ids.size(), MAX_STR_LENGTH);
1626 
1627  std::vector<exII::ex_set> sets(side_boundary_ids.size());
1628 
1629  for (auto i : index_range(side_boundary_ids))
1630  {
1631  boundary_id_type ss_id = side_boundary_ids[i];
1633 
1634  sets[i].id = ss_id;
1635  sets[i].type = exII::EX_SIDE_SET;
1636  sets[i].num_entry = elem[ss_id].size();
1637  sets[i].num_distribution_factor = 0;
1638  sets[i].entry_list = elem[ss_id].data();
1639  sets[i].extra_list = side[ss_id].data();
1640  sets[i].distribution_factor_list = nullptr;
1641  }
1642 
1643  ex_err = exII::ex_put_sets(ex_id, side_boundary_ids.size(), sets.data());
1644  EX_CHECK_ERR(ex_err, "Error writing sidesets");
1645 
1646  ex_err = exII::ex_put_names(ex_id, exII::EX_SIDE_SET, names_table.get_char_star_star());
1647  EX_CHECK_ERR(ex_err, "Error writing sideset names");
1648  }
1649 }
1650 
1651 
1652 
1654 {
1655  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1656  return;
1657 
1658  // Maps from nodeset id to the nodes
1659  std::map<boundary_id_type, std::vector<int>> node;
1660 
1661  // Accumulate the vectors to pass into ex_put_node_set
1662  // build_node_list() builds a list of (node-id, bc-id) tuples.
1663  for (const auto & t : mesh.get_boundary_info().build_node_list())
1664  node[std::get<1>(t)].push_back(std::get<0>(t) + 1);
1665 
1666  std::vector<boundary_id_type> node_boundary_ids;
1667  mesh.get_boundary_info().build_node_boundary_ids(node_boundary_ids);
1668 
1669  // Write out the nodeset names, but only if there is something to write
1670  if (node_boundary_ids.size() > 0)
1671  {
1672  NamesData names_table(node_boundary_ids.size(), MAX_STR_LENGTH);
1673 
1674  for (const auto & nodeset_id : node_boundary_ids)
1675  {
1676  int actual_id = nodeset_id;
1677 
1678  names_table.push_back_entry(mesh.get_boundary_info().get_nodeset_name(nodeset_id));
1679 
1680  ex_err = exII::ex_put_node_set_param(ex_id, actual_id, node[nodeset_id].size(), 0);
1681  EX_CHECK_ERR(ex_err, "Error writing nodeset parameters");
1682 
1683  ex_err = exII::ex_put_node_set(ex_id, actual_id, node[nodeset_id].data());
1684  EX_CHECK_ERR(ex_err, "Error writing nodesets");
1685  }
1686 
1687  // Write out the nodeset names
1688  ex_err = exII::ex_put_names(ex_id, exII::EX_NODE_SET, names_table.get_char_star_star());
1689  EX_CHECK_ERR(ex_err, "Error writing nodeset names");
1690  }
1691 }
1692 
1693 
1694 
1695 void ExodusII_IO_Helper::initialize_element_variables(std::vector<std::string> names,
1696  const std::vector<std::set<subdomain_id_type>> & vars_active_subdomains)
1697 {
1698  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1699  return;
1700 
1701  // Quick return if there are no element variables to write
1702  if (names.size() == 0)
1703  return;
1704 
1705  // Quick return if we have already called this function
1707  return;
1708 
1709  // Be sure that variables in the file match what we are asking for
1710  if (num_elem_vars > 0)
1711  {
1712  this->check_existing_vars(ELEMENTAL, names, this->elem_var_names);
1713  return;
1714  }
1715 
1716  // Set the flag so we can skip this stuff on subsequent calls to
1717  // initialize_element_variables()
1718  _elem_vars_initialized = true;
1719 
1720  this->write_var_names(ELEMENTAL, names);
1721 
1722  // Use the truth table to indicate which subdomain/variable pairs are
1723  // active according to vars_active_subdomains.
1724  std::vector<int> truth_tab(num_elem_blk*num_elem_vars, 0);
1725  for (auto var_num : index_range(vars_active_subdomains))
1726  {
1727  // If the list of active subdomains is empty, it is interpreted as being
1728  // active on *all* subdomains.
1729  std::set<subdomain_id_type> current_set;
1730  if (vars_active_subdomains[var_num].empty())
1731  for (auto block_id : block_ids)
1732  current_set.insert(cast_int<subdomain_id_type>(block_id));
1733  else
1734  current_set = vars_active_subdomains[var_num];
1735 
1736  // Find index into the truth table for each id in current_set.
1737  for (auto block_id : current_set)
1738  {
1739  auto it = std::find(block_ids.begin(), block_ids.end(), block_id);
1740  if (it == block_ids.end())
1741  libmesh_error_msg("ExodusII_IO_Helper: block id " << block_id << " not found in block_ids.");
1742 
1743  std::size_t block_index =
1744  std::distance(block_ids.begin(), it);
1745 
1746  std::size_t truth_tab_index = block_index*num_elem_vars + var_num;
1747  truth_tab[truth_tab_index] = 1;
1748  }
1749  }
1750 
1752  num_elem_blk,
1753  num_elem_vars,
1754  truth_tab.data());
1755  EX_CHECK_ERR(ex_err, "Error writing element truth table.");
1756 }
1757 
1758 
1759 
1760 void ExodusII_IO_Helper::initialize_nodal_variables(std::vector<std::string> names)
1761 {
1762  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1763  return;
1764 
1765  // Quick return if there are no nodal variables to write
1766  if (names.size() == 0)
1767  return;
1768 
1769  // Quick return if we have already called this function
1771  return;
1772 
1773  // Be sure that variables in the file match what we are asking for
1774  if (num_nodal_vars > 0)
1775  {
1776  this->check_existing_vars(NODAL, names, this->nodal_var_names);
1777  return;
1778  }
1779 
1780  // Set the flag so we can skip the rest of this function on subsequent calls.
1781  _nodal_vars_initialized = true;
1782 
1783  this->write_var_names(NODAL, names);
1784 }
1785 
1786 
1787 
1788 void ExodusII_IO_Helper::initialize_global_variables(std::vector<std::string> names)
1789 {
1790  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1791  return;
1792 
1793  // Quick return if there are no global variables to write
1794  if (names.size() == 0)
1795  return;
1796 
1798  return;
1799 
1800  // Be sure that variables in the file match what we are asking for
1801  if (num_global_vars > 0)
1802  {
1803  this->check_existing_vars(GLOBAL, names, this->global_var_names);
1804  return;
1805  }
1806 
1807  _global_vars_initialized = true;
1808 
1809  this->write_var_names(GLOBAL, names);
1810 }
1811 
1812 
1813 
1815  std::vector<std::string> & names,
1816  std::vector<std::string> & names_from_file)
1817 {
1818  // There may already be global variables in the file (for example,
1819  // if we're appending) and in that case, we
1820  // 1.) Cannot initialize them again.
1821  // 2.) Should check to be sure that the global variable names are the same.
1822 
1823  // Fills up names_from_file for us
1824  this->read_var_names(type);
1825 
1826  // Both the names of the global variables and their order must match
1827  if (names_from_file != names)
1828  {
1829  libMesh::err << "Error! The Exodus file already contains the variables:" << std::endl;
1830  for (const auto & name : names_from_file)
1831  libMesh::out << name << std::endl;
1832 
1833  libMesh::err << "And you asked to write:" << std::endl;
1834  for (const auto & name : names)
1835  libMesh::out << name << std::endl;
1836 
1837  libmesh_error_msg("Cannot overwrite existing variables in Exodus II file.");
1838  }
1839 }
1840 
1841 
1842 
1844 {
1845  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1846  return;
1847 
1848  if (_single_precision)
1849  {
1850  float cast_time = float(time);
1851  ex_err = exII::ex_put_time(ex_id, timestep, &cast_time);
1852  }
1853  else
1854  {
1855  ex_err = exII::ex_put_time(ex_id, timestep, &time);
1856  }
1857  EX_CHECK_ERR(ex_err, "Error writing timestep.");
1858 
1860  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
1861 }
1862 
1863 
1864 
1866  const std::vector<Real> & values,
1867  int timestep,
1868  const std::vector<std::set<subdomain_id_type>> & vars_active_subdomains)
1869 {
1870  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1871  return;
1872 
1873  // Loop over the element blocks and write the data one block at a time
1874  std::map<subdomain_id_type, std::vector<unsigned int>> subdomain_map;
1875 
1876  // Ask the file how many element vars it has, store it in the num_elem_vars variable.
1878  EX_CHECK_ERR(ex_err, "Error reading number of elemental variables.");
1879 
1880  // loop through element and map between block and element vector
1881  for (const auto & elem : mesh.active_element_ptr_range())
1882  subdomain_map[elem->subdomain_id()].push_back(elem->id());
1883 
1884  // Use mesh.n_elem() to access into the values vector rather than
1885  // the number of elements the Exodus writer thinks the mesh has,
1886  // which may not include inactive elements.
1888 
1889  // For each variable, create a 'data' array which holds all the elemental variable
1890  // values *for a given block* on this processor, then write that data vector to file
1891  // before moving onto the next block.
1892  libmesh_assert_equal_to(vars_active_subdomains.size(), static_cast<unsigned>(num_elem_vars));
1893  for (unsigned int i=0; i<static_cast<unsigned>(num_elem_vars); ++i)
1894  {
1895  // The size of the subdomain map is the number of blocks.
1896  std::map<subdomain_id_type, std::vector<unsigned int>>::iterator it = subdomain_map.begin();
1897 
1898  const std::set<subdomain_id_type> & active_subdomains = vars_active_subdomains[i];
1899  for (unsigned int j=0; it!=subdomain_map.end(); ++it, ++j)
1900  {
1901  // Skip any variable/subdomain pairs that are inactive.
1902  // Note that if active_subdomains is empty, it is interpreted
1903  // as being active on *all* subdomains.
1904  if (!active_subdomains.empty())
1905  if (active_subdomains.find(it->first) == active_subdomains.end())
1906  {
1907  continue;
1908  }
1909 
1910  const std::vector<unsigned int> & elem_nums = (*it).second;
1911  const unsigned int num_elems_this_block =
1912  cast_int<unsigned int>(elem_nums.size());
1913  std::vector<Real> data(num_elems_this_block);
1914 
1915  for (unsigned int k=0; k<num_elems_this_block; ++k)
1916  data[k] = values[i*n_elem + elem_nums[k]];
1917 
1918  if (_single_precision)
1919  {
1920  std::vector<float> cast_data(data.begin(), data.end());
1921 
1923  timestep,
1924  i+1,
1925  this->get_block_id(j),
1926  num_elems_this_block,
1927  cast_data.data());
1928  }
1929  else
1930  {
1932  timestep,
1933  i+1,
1934  this->get_block_id(j),
1935  num_elems_this_block,
1936  data.data());
1937  }
1938  EX_CHECK_ERR(ex_err, "Error writing element values.");
1939  }
1940  }
1941 
1943  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
1944 }
1945 
1946 
1947 
1948 void ExodusII_IO_Helper::write_nodal_values(int var_id, const std::vector<Real> & values, int timestep)
1949 {
1950  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1951  return;
1952 
1953  if (_single_precision)
1954  {
1955  std::vector<float> cast_values(values.begin(), values.end());
1956  ex_err = exII::ex_put_nodal_var(ex_id, timestep, var_id, num_nodes, cast_values.data());
1957  }
1958  else
1959  {
1960  ex_err = exII::ex_put_nodal_var(ex_id, timestep, var_id, num_nodes, values.data());
1961  }
1962  EX_CHECK_ERR(ex_err, "Error writing nodal values.");
1963 
1965  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
1966 }
1967 
1968 
1969 
1970 void ExodusII_IO_Helper::write_information_records(const std::vector<std::string> & records)
1971 {
1972  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1973  return;
1974 
1975  // There may already be information records in the file (for
1976  // example, if we're appending) and in that case, according to the
1977  // Exodus documentation, writing more information records is not
1978  // supported.
1979  int num_info = inquire(exII::EX_INQ_INFO, "Error retrieving the number of information records from file!");
1980  if (num_info > 0)
1981  {
1982  libMesh::err << "Warning! The Exodus file already contains information records.\n"
1983  << "Exodus does not support writing additional records in this situation."
1984  << std::endl;
1985  return;
1986  }
1987 
1988  int num_records = cast_int<int>(records.size());
1989 
1990  if (num_records > 0)
1991  {
1992  NamesData info(num_records, MAX_LINE_LENGTH);
1993 
1994  // If an entry is longer than MAX_LINE_LENGTH characters it's not an error, we just
1995  // write the first MAX_LINE_LENGTH characters to the file.
1996  for (const auto & record : records)
1997  info.push_back_entry(record);
1998 
1999  ex_err = exII::ex_put_info(ex_id, num_records, info.get_char_star_star());
2000  EX_CHECK_ERR(ex_err, "Error writing global values.");
2001 
2003  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
2004  }
2005 }
2006 
2007 
2008 
2009 void ExodusII_IO_Helper::write_global_values(const std::vector<Real> & values, int timestep)
2010 {
2011  if ((_run_only_on_proc0) && (this->processor_id() != 0))
2012  return;
2013 
2014  if (_single_precision)
2015  {
2016  std::vector<float> cast_values(values.begin(), values.end());
2017  ex_err = exII::ex_put_glob_vars(ex_id, timestep, num_global_vars, cast_values.data());
2018  }
2019  else
2020  {
2021  ex_err = exII::ex_put_glob_vars(ex_id, timestep, num_global_vars, values.data());
2022  }
2023  EX_CHECK_ERR(ex_err, "Error writing global values.");
2024 
2026  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
2027 }
2028 
2029 
2030 
2031 void ExodusII_IO_Helper::read_global_values(std::vector<Real> & values, int timestep)
2032 {
2033  if ((_run_only_on_proc0) && (this->processor_id() != 0))
2034  return;
2035 
2036  values.clear();
2037  values.resize(num_global_vars);
2038  ex_err = exII::ex_get_glob_vars(ex_id, timestep, num_global_vars, values.data());
2039  EX_CHECK_ERR(ex_err, "Error reading global values.");
2040 }
2041 
2042 
2043 
2045 {
2047 }
2048 
2049 
2050 
2052 {
2053  _write_as_dimension = dim;
2054 }
2055 
2056 
2057 
2059 {
2060  _coordinate_offset = p;
2061 }
2062 
2063 
2064 std::vector<std::string> ExodusII_IO_Helper::get_complex_names(const std::vector<std::string> & names) const
2065 {
2066  std::vector<std::string>::const_iterator names_it = names.begin();
2067  std::vector<std::string>::const_iterator names_end = names.end();
2068 
2069  std::vector<std::string> complex_names;
2070 
2071  // This will loop over all names and create new "complex" names
2072  // (i.e. names that start with r_, i_ or a_
2073  for (; names_it != names_end; ++names_it)
2074  {
2075  std::stringstream name_real, name_imag, name_abs;
2076  name_real << "r_" << *names_it;
2077  name_imag << "i_" << *names_it;
2078  name_abs << "a_" << *names_it;
2079 
2080  complex_names.push_back(name_real.str());
2081  complex_names.push_back(name_imag.str());
2082  complex_names.push_back(name_abs.str());
2083  }
2084 
2085  return complex_names;
2086 }
2087 
2088 
2089 
2090 std::vector<std::set<subdomain_id_type>> ExodusII_IO_Helper::get_complex_vars_active_subdomains(
2091  const std::vector<std::set<subdomain_id_type>> & vars_active_subdomains) const
2092 {
2093  std::vector<std::set<subdomain_id_type>> complex_vars_active_subdomains;
2094 
2095  for (auto & s : vars_active_subdomains)
2096  {
2097  // Push back the same data three times to match the tripling of the variables
2098  // for the real, imag, and modulus for the complex-valued solution.
2099  complex_vars_active_subdomains.push_back(s);
2100  complex_vars_active_subdomains.push_back(s);
2101  complex_vars_active_subdomains.push_back(s);
2102  }
2103 
2104  return complex_vars_active_subdomains;
2105 }
2106 
2107 
2108 
2109 // ------------------------------------------------------------
2110 // ExodusII_IO_Helper::Conversion class members
2112 {
2113  init_element_equivalence_map();
2114 
2115  // Do only upper-case comparisons
2116  std::transform(type_str.begin(), type_str.end(), type_str.begin(), ::toupper);
2117 
2118  std::map<std::string, ElemType>::iterator it =
2119  element_equivalence_map.find(type_str);
2120 
2121  if (it != element_equivalence_map.end())
2122  return assign_conversion( it->second );
2123  else
2124  libmesh_error_msg("ERROR! Unrecognized element type_str: " << type_str);
2125 }
2126 
2127 
2128 
2130 {
2131  switch (type)
2132  {
2133  case NODEELEM:
2134  {
2135  const Conversion conv(nodeelem_node_map,
2136  ARRAY_LENGTH(nodeelem_node_map),
2137  nodeelem_node_map, // inverse node map same as forward node map
2138  ARRAY_LENGTH(nodeelem_node_map),
2139  nullptr, // NODELEM doesn't have any edges
2140  0,
2141  nullptr,
2142  0,
2143  NODEELEM, "SPHERE");
2144  return conv;
2145  }
2146 
2147  case EDGE2:
2148  {
2149  const Conversion conv(edge2_node_map,
2150  ARRAY_LENGTH(edge2_node_map),
2151  edge2_node_map, // inverse node map same as forward node map
2152  ARRAY_LENGTH(edge2_node_map),
2153  edge_edge_map,
2154  ARRAY_LENGTH(edge_edge_map),
2155  edge_inverse_edge_map,
2156  ARRAY_LENGTH(edge_inverse_edge_map),
2157  EDGE2, "EDGE2");
2158  return conv;
2159  }
2160  case EDGE3:
2161  {
2162  const Conversion conv(edge3_node_map,
2163  ARRAY_LENGTH(edge3_node_map),
2164  edge3_node_map, // inverse node map same as forward node map
2165  ARRAY_LENGTH(edge3_node_map),
2166  edge_edge_map,
2167  ARRAY_LENGTH(edge_edge_map),
2168  edge_inverse_edge_map,
2169  ARRAY_LENGTH(edge_inverse_edge_map),
2170  EDGE3, "EDGE3");
2171  return conv;
2172  }
2173  case QUAD4:
2174  {
2175  const Conversion conv(quad4_node_map,
2176  ARRAY_LENGTH(quad4_node_map),
2177  quad4_node_map, // inverse node map same as forward node map
2178  ARRAY_LENGTH(quad4_node_map),
2179  quad_edge_map,
2180  ARRAY_LENGTH(quad_edge_map),
2181  quad_inverse_edge_map,
2182  ARRAY_LENGTH(quad_inverse_edge_map),
2183  QUAD4,
2184  "QUAD4");
2185  return conv;
2186  }
2187 
2188  case QUADSHELL4:
2189  {
2190  return Conversion(quad4_node_map,
2191  ARRAY_LENGTH(quad4_node_map), // node mapping is the same as for quad4
2192  quad4_node_map,
2193  ARRAY_LENGTH(quad4_node_map),
2194  quadshell4_edge_map,
2195  ARRAY_LENGTH(quadshell4_edge_map),
2196  quadshell4_inverse_edge_map,
2197  ARRAY_LENGTH(quadshell4_inverse_edge_map),
2198  quadshell4_shellface_map,
2199  ARRAY_LENGTH(quadshell4_shellface_map),
2200  quadshell4_inverse_shellface_map,
2201  ARRAY_LENGTH(quadshell4_inverse_shellface_map),
2202  2, // the side index offset for QUADSHELL4 is 2
2203  QUADSHELL4,
2204  "SHELL4");
2205  }
2206 
2207  case QUAD8:
2208  {
2209  const Conversion conv(quad8_node_map,
2210  ARRAY_LENGTH(quad8_node_map),
2211  quad8_node_map, // inverse node map same as forward node map
2212  ARRAY_LENGTH(quad8_node_map),
2213  quad_edge_map,
2214  ARRAY_LENGTH(quad_edge_map),
2215  quad_inverse_edge_map,
2216  ARRAY_LENGTH(quad_inverse_edge_map),
2217  QUAD8,
2218  "QUAD8");
2219  return conv;
2220  }
2221 
2222  case QUADSHELL8:
2223  {
2224  return Conversion(quad8_node_map,
2225  ARRAY_LENGTH(quad8_node_map), // node mapping is the same as for quad8
2226  quad8_node_map,
2227  ARRAY_LENGTH(quad8_node_map),
2228  quadshell4_edge_map,
2229  ARRAY_LENGTH(quadshell4_edge_map),
2230  quadshell4_inverse_edge_map,
2231  ARRAY_LENGTH(quadshell4_inverse_edge_map),
2232  quadshell4_shellface_map,
2233  ARRAY_LENGTH(quadshell4_shellface_map),
2234  quadshell4_inverse_shellface_map,
2235  ARRAY_LENGTH(quadshell4_inverse_shellface_map),
2236  2, // the side index offset for QUADSHELL8 is 2
2237  QUADSHELL8,
2238  "SHELL8");
2239  }
2240 
2241  case QUAD9:
2242  {
2243  const Conversion conv(quad9_node_map,
2244  ARRAY_LENGTH(quad9_node_map),
2245  quad9_node_map, // inverse node map same as forward node map
2246  ARRAY_LENGTH(quad9_node_map),
2247  quad_edge_map,
2248  ARRAY_LENGTH(quad_edge_map),
2249  quad_inverse_edge_map,
2250  ARRAY_LENGTH(quad_inverse_edge_map),
2251  QUAD9,
2252  "QUAD9");
2253  return conv;
2254  }
2255 
2256  case TRI3:
2257  {
2258  return Conversion(tri3_node_map,
2259  ARRAY_LENGTH(tri3_node_map),
2260  tri3_node_map, // inverse node map same as forward node map
2261  ARRAY_LENGTH(tri3_node_map),
2262  tri_edge_map,
2263  ARRAY_LENGTH(tri_edge_map),
2264  tri_inverse_edge_map,
2265  ARRAY_LENGTH(tri_inverse_edge_map),
2266  TRI3,
2267  "TRI3");
2268  }
2269 
2270  case TRISHELL3:
2271  {
2272  return Conversion(tri3_node_map,
2273  ARRAY_LENGTH(tri3_node_map), // node mapping is the same as for tri3
2274  tri3_node_map,
2275  ARRAY_LENGTH(tri3_node_map),
2276  trishell3_edge_map,
2277  ARRAY_LENGTH(trishell3_edge_map),
2278  trishell3_inverse_edge_map,
2279  ARRAY_LENGTH(trishell3_inverse_edge_map),
2280  trishell3_shellface_map,
2281  ARRAY_LENGTH(trishell3_shellface_map),
2282  trishell3_inverse_shellface_map,
2283  ARRAY_LENGTH(trishell3_inverse_shellface_map),
2284  2, // the side index offset for TRISHELL4 is 2
2285  TRISHELL3,
2286  "TRISHELL3");
2287  }
2288 
2289  case TRI3SUBDIVISION:
2290  {
2291  const Conversion conv(tri3_node_map,
2292  ARRAY_LENGTH(tri3_node_map),
2293  tri3_node_map, // inverse node map same as forward node map
2294  ARRAY_LENGTH(tri3_node_map),
2295  tri_edge_map,
2296  ARRAY_LENGTH(tri_edge_map),
2297  tri_inverse_edge_map,
2298  ARRAY_LENGTH(tri_inverse_edge_map),
2300  "TRI3");
2301  return conv;
2302  }
2303 
2304  case TRI6:
2305  {
2306  const Conversion conv(tri6_node_map,
2307  ARRAY_LENGTH(tri6_node_map),
2308  tri6_node_map, // inverse node map same as forward node map
2309  ARRAY_LENGTH(tri6_node_map),
2310  tri_edge_map,
2311  ARRAY_LENGTH(tri_edge_map),
2312  tri_inverse_edge_map,
2313  ARRAY_LENGTH(tri_inverse_edge_map),
2314  TRI6,
2315  "TRI6");
2316  return conv;
2317  }
2318 
2319  case HEX8:
2320  {
2321  const Conversion conv(hex8_node_map,
2322  ARRAY_LENGTH(hex8_node_map),
2323  hex8_node_map, // inverse node map same as forward node map
2324  ARRAY_LENGTH(hex8_node_map),
2325  hex_face_map,
2326  ARRAY_LENGTH(hex_face_map),
2327  hex_inverse_face_map,
2328  ARRAY_LENGTH(hex_inverse_face_map),
2329  HEX8,
2330  "HEX8");
2331  return conv;
2332  }
2333 
2334  case HEX20:
2335  {
2336  const Conversion conv(hex20_node_map,
2337  ARRAY_LENGTH(hex20_node_map),
2338  hex20_node_map, // inverse node map same as forward node map
2339  ARRAY_LENGTH(hex20_node_map),
2340  hex_face_map,
2341  ARRAY_LENGTH(hex_face_map),
2342  hex_inverse_face_map,
2343  ARRAY_LENGTH(hex_inverse_face_map),
2344  HEX20,
2345  "HEX20");
2346  return conv;
2347  }
2348 
2349  case HEX27:
2350  {
2351  const Conversion conv(hex27_node_map,
2352  ARRAY_LENGTH(hex27_node_map),
2353  hex27_inverse_node_map, // different inverse node map for Hex27!
2354  ARRAY_LENGTH(hex27_inverse_node_map),
2355  hex27_face_map,
2356  ARRAY_LENGTH(hex27_face_map),
2357  hex27_inverse_face_map,
2358  ARRAY_LENGTH(hex27_inverse_face_map),
2359  HEX27,
2360  "HEX27");
2361  return conv;
2362  }
2363 
2364  case TET4:
2365  {
2366  const Conversion conv(tet4_node_map,
2367  ARRAY_LENGTH(tet4_node_map),
2368  tet4_node_map, // inverse node map same as forward node map
2369  ARRAY_LENGTH(tet4_node_map),
2370  tet_face_map,
2371  ARRAY_LENGTH(tet_face_map),
2372  tet_inverse_face_map,
2373  ARRAY_LENGTH(tet_inverse_face_map),
2374  TET4,
2375  "TETRA4");
2376  return conv;
2377  }
2378 
2379  case TET10:
2380  {
2381  const Conversion conv(tet10_node_map,
2382  ARRAY_LENGTH(tet10_node_map),
2383  tet10_node_map, // inverse node map same as forward node map
2384  ARRAY_LENGTH(tet10_node_map),
2385  tet_face_map,
2386  ARRAY_LENGTH(tet_face_map),
2387  tet_inverse_face_map,
2388  ARRAY_LENGTH(tet_inverse_face_map),
2389  TET10,
2390  "TETRA10");
2391  return conv;
2392  }
2393 
2394  case PRISM6:
2395  {
2396  const Conversion conv(prism6_node_map,
2397  ARRAY_LENGTH(prism6_node_map),
2398  prism6_node_map, // inverse node map same as forward node map
2399  ARRAY_LENGTH(prism6_node_map),
2400  prism_face_map,
2401  ARRAY_LENGTH(prism_face_map),
2402  prism_inverse_face_map,
2403  ARRAY_LENGTH(prism_inverse_face_map),
2404  PRISM6,
2405  "WEDGE");
2406  return conv;
2407  }
2408 
2409  case PRISM15:
2410  {
2411  const Conversion conv(prism15_node_map,
2412  ARRAY_LENGTH(prism15_node_map),
2413  prism15_node_map, // inverse node map same as forward node map
2414  ARRAY_LENGTH(prism15_node_map),
2415  prism_face_map,
2416  ARRAY_LENGTH(prism_face_map),
2417  prism_inverse_face_map,
2418  ARRAY_LENGTH(prism_inverse_face_map),
2419  PRISM15,
2420  "WEDGE15");
2421  return conv;
2422  }
2423 
2424  case PRISM18:
2425  {
2426  const Conversion conv(prism18_node_map,
2427  ARRAY_LENGTH(prism18_node_map),
2428  prism18_node_map, // inverse node map same as forward node map
2429  ARRAY_LENGTH(prism18_node_map),
2430  prism_face_map,
2431  ARRAY_LENGTH(prism_face_map),
2432  prism_inverse_face_map,
2433  ARRAY_LENGTH(prism_inverse_face_map),
2434  PRISM18,
2435  "WEDGE18");
2436  return conv;
2437  }
2438 
2439  case PYRAMID5:
2440  {
2441  const Conversion conv(pyramid5_node_map,
2442  ARRAY_LENGTH(pyramid5_node_map),
2443  pyramid5_node_map, // inverse node map same as forward node map
2444  ARRAY_LENGTH(pyramid5_node_map),
2445  pyramid_face_map,
2446  ARRAY_LENGTH(pyramid_face_map),
2447  pyramid_inverse_face_map,
2448  ARRAY_LENGTH(pyramid_inverse_face_map),
2449  PYRAMID5,
2450  "PYRAMID5");
2451  return conv;
2452  }
2453 
2454  case PYRAMID13:
2455  {
2456  const Conversion conv(pyramid13_node_map,
2457  ARRAY_LENGTH(pyramid13_node_map),
2458  pyramid13_node_map, // inverse node map same as forward node map
2459  ARRAY_LENGTH(pyramid13_node_map),
2460  pyramid_face_map,
2461  ARRAY_LENGTH(pyramid_face_map),
2462  pyramid_inverse_face_map,
2463  ARRAY_LENGTH(pyramid_inverse_face_map),
2464  PYRAMID13,
2465  "PYRAMID13");
2466  return conv;
2467  }
2468 
2469  case PYRAMID14:
2470  {
2471  const Conversion conv(pyramid14_node_map,
2472  ARRAY_LENGTH(pyramid14_node_map),
2473  pyramid14_node_map, // inverse node map same as forward node map
2474  ARRAY_LENGTH(pyramid14_node_map),
2475  pyramid_face_map,
2476  ARRAY_LENGTH(pyramid_face_map),
2477  pyramid_inverse_face_map,
2478  ARRAY_LENGTH(pyramid_inverse_face_map),
2479  PYRAMID14,
2480  "PYRAMID14");
2481  return conv;
2482  }
2483 
2484  default:
2485  libmesh_error_msg("Unsupported element type: " << type);
2486  }
2487 }
2488 
2489 
2490 
2492 {
2493  // If we asked for a side that doesn't exist, return an invalid_id
2494  // and allow higher-level code to handle it.
2495  if (static_cast<size_t>(i) >= side_map_size)
2496  return invalid_id;
2497 
2498  return side_map[i];
2499 }
2500 
2501 
2502 
2503 ExodusII_IO_Helper::NamesData::NamesData(size_t n_strings, size_t string_length) :
2504  data_table(n_strings),
2505  data_table_pointers(n_strings),
2506  counter(0),
2507  table_size(n_strings)
2508 {
2509  for (size_t i=0; i<n_strings; ++i)
2510  {
2511  data_table[i].resize(string_length + 1);
2512 
2513  // Properly terminate these C-style strings, just to be safe.
2514  data_table[i][0] = '\0';
2515 
2516  // Set pointer into the data_table
2517  data_table_pointers[i] = data_table[i].data();
2518  }
2519 }
2520 
2521 
2522 
2524 {
2525  libmesh_assert_less (counter, table_size);
2526 
2527  // 1.) Copy the C++ string into the vector<char>...
2528  size_t num_copied = name.copy(data_table[counter].data(), data_table[counter].size()-1);
2529 
2530  // 2.) ...And null-terminate it.
2531  data_table[counter][num_copied] = '\0';
2532 
2533  // Go to next row
2534  ++counter;
2535 }
2536 
2537 
2538 
2540 {
2541  return data_table_pointers.data();
2542 }
2543 
2544 
2545 
2547 {
2548  if (static_cast<unsigned>(i) >= table_size)
2549  libmesh_error_msg("Requested char * " << i << " but only have " << table_size << "!");
2550 
2551  else
2552  return data_table[i].data();
2553 }
2554 
2555 
2556 } // namespace libMesh
2557 
2558 
2559 
2560 #endif // #ifdef LIBMESH_HAVE_EXODUS_API
std::string name(const ElemQuality q)
Definition: elem_quality.C:42
virtual void write_sidesets(const MeshBase &mesh)
void use_mesh_dimension_instead_of_spatial_dimension(bool val)
EXODUS_EXPORT int ex_put_glob_vars(int exoid, int time_step, int num_glob_vars, const void *glob_var_vals)
EXODUS_EXPORT int ex_get_elem_block(int exoid, ex_entity_id elem_blk_id, char *elem_type, void_int *num_elem_this_blk, void_int *num_nodes_per_elem, void_int *num_attr)
EXODUS_EXPORT int ex_get_side_set_ids(int exoid, void_int *ids)
EXODUS_EXPORT int ex_put_info(int exoid, int num_info, char *info[])
std::vector< Real > nodal_var_values
std::vector< std::string > get_complex_names(const std::vector< std::string > &names) const
EXODUS_EXPORT int ex_get_name(int exoid, ex_entity_type obj_type, ex_entity_id entity_id, char *name)
std::vector< int > num_sides_per_set
const char * get_elem_type() const
std::vector< std::string > elem_var_names
virtual dof_id_type n_active_elem() const =0
A geometric point in (x,y,z) space associated with a DOF.
Definition: node.h:52
EXODUS_EXPORT int ex_get_qa(int exoid, char *qa_record[][4])
EXODUS_EXPORT int ex_get_node_set_param(int exoid, ex_entity_id node_set_id, void_int *num_nodes_in_set, void_int *num_df_in_set)
void write_var_names_impl(const char *var_type, int &count, std::vector< std::string > &names)
EXODUS_EXPORT int ex_get_var_names(int exoid, const char *var_type, int num_vars, char *var_names[])
EXODUS_EXPORT int ex_get_node_num_map(int exoid, void_int *node_map)
EXODUS_EXPORT int ex_get_side_set(int exoid, ex_entity_id side_set_id, void_int *side_set_elem_list, void_int *side_set_side_list)
void write_var_names(ExodusVarType type, std::vector< std::string > &names)
std::map< int, std::string > id_to_ss_names
std::string get_block_name(int index)
void read_nodal_var_values(std::string nodal_var_name, int time_step)
ExodusII_IO_Helper(const ParallelObject &parent, bool v=false, bool run_only_on_proc0=true, bool single_precision=false)
dof_id_type n_elem(const MeshBase::const_element_iterator &begin, const MeshBase::const_element_iterator &end)
Definition: mesh_tools.C:702
void write_information_records(const std::vector< std::string > &records)
void read_var_names_impl(const char *var_type, int &count, std::vector< std::string > &result)
void write_as_dimension(unsigned dim)
EXODUS_EXPORT int ex_put_concat_elem_block(int exoid, const void_int *elem_blk_id, char *elem_type[], const void_int *num_elem_this_blk, const void_int *num_nodes_per_elem, const void_int *num_attr, int define_maps)
void print_nodes(std::ostream &out=libMesh::out)
std::string get_side_set_name(int index)
EXODUS_EXPORT int ex_get_node_set(int exoid, ex_entity_id node_set_id, void_int *node_set_node_list)
unsigned short int side
Definition: xdr_io.C:50
EXODUS_EXPORT int ex_put_node_set_param(int exoid, ex_entity_id node_set_id, int64_t num_nodes_in_set, int64_t num_dist_in_set)
The base class for all geometric element types.
Definition: elem.h:100
MeshBase & mesh
IntRange< std::size_t > index_range(const std::vector< T > &vec)
Definition: int_range.h:104
ExodusII_IO_Helper::Conversion assign_conversion(std::string type_str)
EXODUS_EXPORT int ex_get_elem_num_map(int exoid, void_int *elem_map)
int inquire(int req_info, std::string error_msg="")
virtual void write_nodesets(const MeshBase &mesh)
EXODUS_EXPORT int ex_put_var_param(int exoid, const char *var_type, int num_vars)
EXODUS_EXPORT int ex_put_coord(int exoid, const void *x_coor, const void *y_coor, const void *z_coor)
void build_side_boundary_ids(std::vector< boundary_id_type > &b_ids) const
virtual SimpleRange< element_iterator > active_element_ptr_range()=0
std::vector< Real > time_steps
virtual void write_elements(const MeshBase &mesh, bool use_discontinuous=false)
std::vector< std::vector< char > > data_table
EXODUS_EXPORT int ex_get_node_set_ids(int exoid, void_int *ids)
const BoundaryInfo & get_boundary_info() const
Definition: mesh_base.h:131
long double max(long double a, double b)
void write_nodal_values(int var_id, const std::vector< Real > &values, int timestep)
void check_existing_vars(ExodusVarType type, std::vector< std::string > &names, std::vector< std::string > &names_from_file)
EXODUS_EXPORT int ex_put_init(int exoid, const char *title, int64_t num_dim, int64_t num_nodes, int64_t num_elem, int64_t num_elem_blk, int64_t num_node_sets, int64_t num_side_sets)
EXODUS_EXPORT int ex_put_var_names(int exoid, const char *var_type, int num_vars, char *var_names[])
Base class for Mesh.
Definition: mesh_base.h:77
void write_element_values(const MeshBase &mesh, const std::vector< Real > &values, int timestep, const std::vector< std::set< subdomain_id_type >> &vars_active_subdomains)
void write_timestep(int timestep, Real time)
std::map< int, int > libmesh_elem_num_to_exodus
EXODUS_EXPORT int ex_put_sets(int exoid, size_t set_count, const struct ex_set *sets)
void build_side_list(std::vector< dof_id_type > &element_id_list, std::vector< unsigned short int > &side_list, std::vector< boundary_id_type > &bc_id_list) const
EXODUS_EXPORT int ex_inquire(int exoid, int inquiry, void_int *, float *, char *)
EXODUS_EXPORT int ex_get_init(int exoid, char *title, void_int *num_dim, void_int *num_nodes, void_int *num_elem, void_int *num_elem_blk, void_int *num_node_sets, void_int *num_side_sets)
EXODUS_EXPORT int ex_get_all_times(int exoid, void *time_values)
EXODUS_EXPORT int ex_put_elem_var(int exoid, int time_step, int elem_var_index, ex_entity_id elem_blk_id, int64_t num_elem_this_blk, const void *elem_var_vals)
virtual node_iterator nodes_begin()=0
EXODUS_EXPORT int ex_put_node_set(int exoid, ex_entity_id node_set_id, const void_int *node_set_node_list)
virtual void initialize_element_variables(std::vector< std::string > names, const std::vector< std::set< subdomain_id_type >> &vars_active_subdomains)
void open(const char *filename, bool read_only)
void build_node_list(std::vector< dof_id_type > &node_id_list, std::vector< boundary_id_type > &bc_id_list) const
int8_t boundary_id_type
Definition: id_types.h:51
void active_family_tree_by_side(std::vector< const Elem *> &family, const unsigned int side, const bool reset=true) const
Definition: elem.C:1630
void push_back_entry(const std::string &name)
dof_id_type id() const
Definition: dof_object.h:655
void read_global_values(std::vector< Real > &values, int timestep)
virtual unsigned int n_nodes() const =0
EXODUS_EXPORT int ex_get_elem_blk_ids(int exoid, void_int *ids)
std::vector< std::string > global_var_names
EXODUS_EXPORT int ex_get_nodal_var(int exoid, int time_step, int nodal_var_index, int64_t num_nodes, void *nodal_var_vals)
virtual SimpleRange< node_iterator > node_ptr_range()=0
void read_var_names(ExodusVarType type)
EXODUS_EXPORT int ex_get_side_set_param(int exoid, ex_entity_id side_set_id, void_int *num_side_in_set, void_int *num_dist_fact_in_set)
std::string & subdomain_name(subdomain_id_type id)
Definition: mesh_base.C:538
OStreamProxy err(std::cerr)
void write_global_values(const std::vector< Real > &values, int timestep)
EXODUS_EXPORT int ex_put_names(int exoid, ex_entity_type obj_type, char *names[])
void build_shellface_list(std::vector< dof_id_type > &element_id_list, std::vector< unsigned short int > &shellface_list, std::vector< boundary_id_type > &bc_id_list) const
An object whose state is distributed along a set of processors.
const std::string & get_nodeset_name(boundary_id_type id) const
std::map< int, int > libmesh_node_num_to_exodus
void initialize_global_variables(std::vector< std::string > names)
EXODUS_EXPORT int ex_get_elem_var(int exoid, int time_step, int elem_var_index, ex_entity_id elem_blk_id, int64_t num_elem_this_blk, void *elem_var_vals)
std::string enum_to_string(const T e)
std::vector< int > num_node_df_per_set
virtual const Elem * elem_ptr(const dof_id_type i) const =0
void build_node_boundary_ids(std::vector< boundary_id_type > &b_ids) const
EXODUS_EXPORT int ex_get_coord(int exoid, void *x_coor, void *y_coor, void *z_coor)
void message(const std::string &msg)
virtual void initialize(std::string title, const MeshBase &mesh, bool use_discontinuous=false)
EXODUS_EXPORT int ex_get_glob_vars(int exoid, int time_step, int num_glob_vars, void *glob_var_vals)
virtual void create(std::string filename)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual node_iterator nodes_end()=0
EXODUS_EXPORT int ex_put_node_num_map(int exoid, const void_int *node_map)
EXODUS_EXPORT int ex_close(int exoid)
unsigned int spatial_dimension() const
Definition: mesh_base.C:135
void read_sideset(int id, int offset)
EXODUS_EXPORT int ex_put_elem_conn(int exoid, ex_entity_id elem_blk_id, const void_int *connect)
NamesData(size_t n_strings, size_t string_length)
const std::string & get_sideset_name(boundary_id_type id) const
virtual const Elem & elem_ref(const dof_id_type i) const
Definition: mesh_base.h:504
std::map< int, std::string > id_to_block_names
EXODUS_EXPORT int ex_get_elem_conn(int exoid, ex_entity_id elem_blk_id, void_int *connect)
unsigned int mesh_dimension() const
Definition: mesh_base.C:126
EXODUS_EXPORT int ex_get_var_param(int exoid, const char *var_type, int *num_vars)
std::map< int, std::string > id_to_ns_names
void initialize_nodal_variables(std::vector< std::string > names)
std::vector< int > num_nodes_per_set
EXODUS_EXPORT int ex_put_elem_var_tab(int exoid, int num_elem_blk, int num_elem_var, int *elem_var_tab)
std::vector< std::string > nodal_var_names
EXODUS_EXPORT int ex_update(int exoid)
EXODUS_EXPORT void ex_err(const char *module_name, const char *message, int err_num)
std::string get_node_set_name(int index)
IterBase * data
virtual dof_id_type n_elem() const =0
processor_id_type processor_id() const
EXODUS_EXPORT int ex_put_elem_num_map(int exoid, const void_int *elem_map)
OStreamProxy out(std::cout)
void read_elemental_var_values(std::string elemental_var_name, int time_step, std::map< dof_id_type, Real > &elem_var_value_map)
EXODUS_EXPORT int ex_put_time(int exoid, int time_step, const void *time_value)
virtual ElemType type() const =0
long double min(long double a, double b)
A geometric point in (x,y,z) space.
Definition: point.h:38
dof_id_type node_id(const unsigned int i) const
Definition: elem.h:1914
std::vector< int > num_df_per_set
void build_shellface_boundary_ids(std::vector< boundary_id_type > &b_ids) const
virtual void write_nodal_coordinates(const MeshBase &mesh, bool use_discontinuous=false)
uint8_t dof_id_type
Definition: id_types.h:64
EXODUS_EXPORT int ex_put_nodal_var(int exoid, int time_step, int nodal_var_index, int64_t num_nodes, const void *nodal_var_vals)
std::vector< std::set< subdomain_id_type > > get_complex_vars_active_subdomains(const std::vector< std::set< subdomain_id_type >> &vars_active_subdomains) const