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