libMesh::ExodusII_IO_Helper Class Reference

#include <exodusII_io_helper.h>

Inheritance diagram for libMesh::ExodusII_IO_Helper:

Classes

class  Conversion
 
class  ElementMaps
 
class  NamesData
 

Public Member Functions

 ExodusII_IO_Helper (const ParallelObject &parent, bool v=false, bool run_only_on_proc0=true)
 
virtual ~ExodusII_IO_Helper ()
 
const char * get_elem_type () const
 
void open (const char *filename, bool read_only)
 
void read_header ()
 
void print_header ()
 
void read_nodes ()
 
void read_node_num_map ()
 
void print_nodes (std::ostream &out=libMesh::out)
 
void read_block_info ()
 
int get_block_id (int index)
 
std::string get_block_name (int index)
 
int get_side_set_id (int index)
 
std::string get_side_set_name (int index)
 
int get_node_set_id (int index)
 
std::string get_node_set_name (int index)
 
void read_elem_in_block (int block)
 
void read_elem_num_map ()
 
void read_sideset_info ()
 
void read_nodeset_info ()
 
void read_sideset (int id, int offset)
 
void read_nodeset (int id)
 
void close ()
 
int inquire (int req_info, std::string error_msg="")
 
void read_time_steps ()
 
void read_num_time_steps ()
 
void read_nodal_var_values (std::string nodal_var_name, int time_step)
 
void read_elemental_var_values (std::string elemental_var_name, int time_step)
 
virtual void create (std::string filename)
 
virtual void initialize (std::string title, const MeshBase &mesh)
 
void initialize_discontinuous (std::string title, const MeshBase &mesh)
 
virtual void write_nodal_coordinates (const MeshBase &mesh)
 
void write_nodal_coordinates_discontinuous (const MeshBase &mesh)
 
virtual void write_elements (const MeshBase &mesh)
 
void write_elements_discontinuous (const MeshBase &mesh)
 
virtual void write_sidesets (const MeshBase &mesh)
 
virtual void write_nodesets (const MeshBase &mesh)
 
void initialize_element_variables (std::vector< std::string > names)
 
void initialize_nodal_variables (std::vector< std::string > names)
 
void initialize_global_variables (std::vector< std::string > names)
 
void write_timestep (int timestep, Real time)
 
void write_element_values (const MeshBase &mesh, const std::vector< Number > &values, int timestep)
 
void write_nodal_values (int var_id, const std::vector< Number > &values, int timestep)
 
void write_information_records (const std::vector< std::string > &records)
 
void write_global_values (const std::vector< Number > &values, int timestep)
 
void use_mesh_dimension_instead_of_spatial_dimension (bool val)
 
void set_coordinate_offset (Point p)
 
void message (const std::string msg)
 
void message (const std::string msg, int i)
 
const Parallel::Communicatorcomm () const
 
processor_id_type n_processors () const
 
processor_id_type processor_id () const
 

Public Attributes

int ex_id
 
int ex_err
 
int num_dim
 
int num_global_vars
 
int num_nodes
 
int num_elem
 
int num_elem_blk
 
int num_node_sets
 
int num_side_sets
 
int num_elem_this_blk
 
int num_nodes_per_elem
 
int num_attr
 
int num_elem_all_sidesets
 
std::vector< int > block_ids
 
std::vector< int > connect
 
std::vector< int > ss_ids
 
std::vector< int > nodeset_ids
 
std::vector< int > num_sides_per_set
 
std::vector< int > num_nodes_per_set
 
std::vector< int > num_df_per_set
 
std::vector< int > num_node_df_per_set
 
std::vector< int > elem_list
 
std::vector< int > side_list
 
std::vector< int > node_list
 
std::vector< int > id_list
 
std::vector< int > node_num_map
 
std::vector< int > elem_num_map
 
std::vector< Realx
 
std::vector< Realy
 
std::vector< Realz
 
std::vector< char > title
 
std::vector< char > elem_type
 
std::map< int, int > libmesh_elem_num_to_exodus
 
std::vector< int > exodus_elem_num_to_libmesh
 
std::map< int, int > libmesh_node_num_to_exodus
 
std::vector< int > exodus_node_num_to_libmesh
 
int num_time_steps
 
std::vector< Realtime_steps
 
int num_nodal_vars
 
std::vector< std::string > nodal_var_names
 
std::vector< Realnodal_var_values
 
int num_elem_vars
 
std::vector< std::string > elem_var_names
 
std::vector< Realelem_var_values
 
std::vector< std::string > global_var_names
 
std::map< int, std::string > id_to_block_names
 
std::map< int, std::string > id_to_ss_names
 
std::map< int, std::string > id_to_ns_names
 
bool verbose
 
bool opened_for_writing
 
bool opened_for_reading
 
std::string current_filename
 

Protected Attributes

bool _run_only_on_proc0
 
bool _elem_vars_initialized
 
bool _global_vars_initialized
 
bool _nodal_vars_initialized
 
bool _use_mesh_dimension_instead_of_spatial_dimension
 
Point _coordinate_offset
 
const Parallel::Communicator_communicator
 

Private Types

enum  ExodusVarType { NODAL =0, ELEMENTAL =1, GLOBAL =2 }
 

Private Member Functions

void read_var_names (ExodusVarType type)
 
void write_var_names (ExodusVarType type, std::vector< std::string > &names)
 
void check_existing_vars (ExodusVarType type, std::vector< std::string > &names, std::vector< std::string > &names_from_file)
 
void read_var_names_impl (const char *var_type, int &count, std::vector< std::string > &result)
 
void write_var_names_impl (const char *var_type, int &count, std::vector< std::string > &names)
 

Detailed Description

This is the ExodusII_IO_Helper class. This class hides the implementation details of interfacing with the Exodus binary format.

Author
Johw W. Peterson, 2002.

Definition at line 61 of file exodusII_io_helper.h.

Member Enumeration Documentation

Wraps calls to exII::ex_get_var_names() and exII::ex_get_var_param(). The enumeration controls whether nodal, elemental, or global variable names are read and which class members are filled in. NODAL: num_nodal_vars nodal_var_names ELEMENTAL: num_elem_vars elem_var_names GLOBAL: num_global_vars global_var_names

Enumerator
NODAL 
ELEMENTAL 
GLOBAL 

Definition at line 557 of file exodusII_io_helper.h.

557 {NODAL=0, ELEMENTAL=1, GLOBAL=2};

Constructor & Destructor Documentation

libMesh::ExodusII_IO_Helper::ExodusII_IO_Helper ( const ParallelObject parent,
bool  v = false,
bool  run_only_on_proc0 = true 
)

Constructor. Automatically initializes all the private members of the class. Also allows you to set the verbosity level to v=true (on) or v=false (off). The second argument, if true, tells the class to only perform its actions if running on processor zero. If you initialize this to false, the writing methods will run on all processors instead.

Definition at line 225 of file exodusII_io_helper.C.

References elem_type, and title.

227  :
228  ParallelObject(parent),
229  ex_id(0),
230  ex_err(0),
231  num_dim(0),
232  num_global_vars(0),
233  num_nodes(0),
234  num_elem(0),
235  num_elem_blk(0),
236  num_node_sets(0),
237  num_side_sets(0),
240  num_attr(0),
242  num_time_steps(0),
243  num_nodal_vars(0),
244  num_elem_vars(0),
245  verbose(v),
246  opened_for_writing(false),
247  opened_for_reading(false),
248  _run_only_on_proc0(run_only_on_proc0),
249  _elem_vars_initialized(false),
253  {
254  title.resize(MAX_LINE_LENGTH+1);
255  elem_type.resize(MAX_STR_LENGTH);
256  }
libMesh::ExodusII_IO_Helper::~ExodusII_IO_Helper ( )
virtual

Destructor.

Definition at line 260 of file exodusII_io_helper.C.

261 {
262 }

Member Function Documentation

void libMesh::ExodusII_IO_Helper::check_existing_vars ( ExodusVarType  type,
std::vector< std::string > &  names,
std::vector< std::string > &  names_from_file 
)
private

When appending: during initialization, check that variable names in the file match those you attempt to initialize with.

Definition at line 1654 of file exodusII_io_helper.C.

References libMesh::err, libMesh::out, and read_var_names().

Referenced by initialize_element_variables(), initialize_global_variables(), and initialize_nodal_variables().

1657 {
1658  // There may already be global variables in the file (for example,
1659  // if we're appending) and in that case, we
1660  // 1.) Cannot initialize them again.
1661  // 2.) Should check to be sure that the global variable names are the same.
1662 
1663  // Fills up names_from_file for us
1664  this->read_var_names(type);
1665 
1666  // Both the names of the global variables and their order must match
1667  if (names_from_file != names)
1668  {
1669  libMesh::err << "Error! The Exodus file already contains the variables:" << std::endl;
1670  for (unsigned i=0; i<names_from_file.size(); ++i)
1671  libMesh::out << names_from_file[i] << std::endl;
1672 
1673  libMesh::err << "And you asked to write:" << std::endl;
1674  for (unsigned i=0; i<names.size(); ++i)
1675  libMesh::out << names[i] << std::endl;
1676 
1677  libmesh_error();
1678  }
1679 }
void libMesh::ExodusII_IO_Helper::close ( )

Closes the ExodusII mesh file.

Definition at line 688 of file exodusII_io_helper.C.

References _run_only_on_proc0, ex_close(), ex_err, ex_id, message(), and libMesh::ParallelObject::processor_id().

Referenced by libMesh::ExodusII_IO::~ExodusII_IO(), and libMesh::Nemesis_IO_Helper::~Nemesis_IO_Helper().

689 {
690  // Always call close on processor 0.
691  // If we're running on multiple processors, i.e. as one of several Nemesis files,
692  // we call close on all processors...
693  if ((this->processor_id() == 0) || (!_run_only_on_proc0))
694  {
696  EX_CHECK_ERR(ex_err, "Error closing Exodus file.");
697  message("Exodus file closed successfully.");
698  }
699 }
const Parallel::Communicator& libMesh::ParallelObject::comm ( ) const
inlineinherited
Returns
a reference to the Parallel::Communicator object used by this mesh.

Definition at line 86 of file parallel_object.h.

References libMesh::ParallelObject::_communicator.

Referenced by libMesh::__libmesh_petsc_diff_solver_jacobian(), libMesh::__libmesh_petsc_diff_solver_monitor(), libMesh::__libmesh_petsc_diff_solver_residual(), libMesh::__libmesh_petsc_snes_jacobian(), libMesh::__libmesh_petsc_snes_residual(), libMesh::MeshRefinement::_coarsen_elements(), libMesh::ExactSolution::_compute_error(), libMesh::MetisPartitioner::_do_partition(), libMesh::ParmetisPartitioner::_do_repartition(), libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::SlepcEigenSolver< T >::_petsc_shell_matrix_get_diagonal(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_get_diagonal(), libMesh::SlepcEigenSolver< T >::_petsc_shell_matrix_mult(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_mult(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_mult_add(), libMesh::EquationSystems::_read_impl(), libMesh::MeshRefinement::_refine_elements(), libMesh::ParallelMesh::add_elem(), libMesh::ImplicitSystem::add_matrix(), libMesh::ParallelMesh::add_node(), libMesh::System::add_vector(), libMesh::UnstructuredMesh::all_second_order(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::FEMSystem::assemble_qoi(), libMesh::MeshCommunication::assign_global_indices(), libMesh::ParmetisPartitioner::assign_partitioning(), libMesh::DofMap::attach_matrix(), libMesh::MeshTools::bounding_box(), libMesh::System::calculate_norm(), libMesh::MeshRefinement::coarsen_elements(), libMesh::Nemesis_IO_Helper::compute_num_global_elem_blocks(), libMesh::Nemesis_IO_Helper::compute_num_global_nodesets(), libMesh::Nemesis_IO_Helper::compute_num_global_sidesets(), libMesh::Problem_Interface::computeF(), libMesh::Problem_Interface::computeJacobian(), libMesh::Problem_Interface::computePreconditioner(), libMesh::MeshTools::correct_node_proc_ids(), libMesh::MeshCommunication::delete_remote_elements(), libMesh::DofMap::distribute_dofs(), DMlibMeshFunction(), DMlibMeshJacobian(), DMLibMeshSetSystem(), DMVariableBounds_libMesh(), libMesh::MeshRefinement::eliminate_unrefined_patches(), libMesh::WeightedPatchRecoveryErrorEstimator::estimate_error(), libMesh::PatchRecoveryErrorEstimator::estimate_error(), libMesh::JumpErrorEstimator::estimate_error(), libMesh::AdjointRefinementEstimator::estimate_error(), libMesh::MeshRefinement::flag_elements_by_elem_fraction(), libMesh::MeshRefinement::flag_elements_by_error_fraction(), libMesh::MeshRefinement::flag_elements_by_nelem_target(), libMesh::for(), libMesh::CondensedEigenSystem::get_eigenpair(), libMesh::ImplicitSystem::get_linear_solver(), libMesh::LocationMap< T >::init(), libMesh::PetscDiffSolver::init(), libMesh::TimeSolver::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::EigenSystem::init_data(), libMesh::EigenSystem::init_matrices(), libMesh::ParmetisPartitioner::initialize(), libMesh::MeshTools::libmesh_assert_valid_dof_ids(), libMesh::ParallelMesh::libmesh_assert_valid_parallel_flags(), libMesh::MeshTools::libmesh_assert_valid_procids< Elem >(), libMesh::MeshTools::libmesh_assert_valid_procids< Node >(), libMesh::MeshTools::libmesh_assert_valid_refinement_flags(), libMesh::MeshRefinement::limit_level_mismatch_at_edge(), libMesh::MeshRefinement::limit_level_mismatch_at_node(), libMesh::MeshRefinement::make_coarsening_compatible(), libMesh::MeshCommunication::make_elems_parallel_consistent(), libMesh::MeshRefinement::make_flags_parallel_consistent(), libMesh::MeshCommunication::make_node_ids_parallel_consistent(), libMesh::MeshCommunication::make_node_proc_ids_parallel_consistent(), libMesh::MeshCommunication::make_nodes_parallel_consistent(), libMesh::MeshRefinement::make_refinement_compatible(), libMesh::FEMSystem::mesh_position_set(), libMesh::MeshSerializer::MeshSerializer(), libMesh::ParallelMesh::n_active_elem(), libMesh::MeshTools::n_active_levels(), libMesh::BoundaryInfo::n_boundary_conds(), libMesh::BoundaryInfo::n_edge_conds(), libMesh::CondensedEigenSystem::n_global_non_condensed_dofs(), libMesh::MeshTools::n_levels(), libMesh::BoundaryInfo::n_nodeset_conds(), libMesh::MeshTools::n_p_levels(), libMesh::ParallelMesh::parallel_max_elem_id(), libMesh::ParallelMesh::parallel_max_node_id(), libMesh::ParallelMesh::parallel_n_elem(), libMesh::ParallelMesh::parallel_n_nodes(), libMesh::Partitioner::partition(), libMesh::Partitioner::partition_unpartitioned_elements(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::MeshBase::prepare_for_use(), libMesh::System::project_vector(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), libMesh::System::read_header(), libMesh::System::read_legacy_data(), libMesh::System::read_SCALAR_dofs(), libMesh::XdrIO::read_serialized_bc_names(), libMesh::XdrIO::read_serialized_bcs(), libMesh::System::read_serialized_blocked_dof_objects(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::System::read_serialized_vector(), libMesh::MeshBase::recalculate_n_partitions(), libMesh::MeshRefinement::refine_and_coarsen_elements(), libMesh::MeshRefinement::refine_elements(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::MeshBase::subdomain_ids(), libMesh::BoundaryInfo::sync(), libMesh::Parallel::sync_element_data_by_parent_id(), libMesh::MeshRefinement::test_level_one(), libMesh::MeshRefinement::test_unflagged(), libMesh::MeshTools::total_weight(), libMesh::CheckpointIO::write(), libMesh::XdrIO::write(), libMesh::UnstructuredMesh::write(), libMesh::LegacyXdrIO::write_mesh(), libMesh::System::write_SCALAR_dofs(), libMesh::XdrIO::write_serialized_bcs(), libMesh::System::write_serialized_blocked_dof_objects(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), libMesh::XdrIO::write_serialized_nodesets(), and libMesh::DivaIO::write_stream().

87  { return _communicator; }
void libMesh::ExodusII_IO_Helper::create ( std::string  filename)
virtual

Opens an ExodusII mesh file named filename for writing.

Reimplemented in libMesh::Nemesis_IO_Helper.

Definition at line 963 of file exodusII_io_helper.C.

References _run_only_on_proc0, current_filename, ex_id, std::min(), opened_for_writing, libMesh::out, libMesh::ParallelObject::processor_id(), libMesh::Real, and verbose.

Referenced by libMesh::ExodusII_IO::write(), and libMesh::ExodusII_IO::write_nodal_data_common().

964 {
965  // If we're processor 0, always create the file.
966  // If we running on all procs, e.g. as one of several Nemesis files, also
967  // call create there.
968  if ((this->processor_id() == 0) || (!_run_only_on_proc0))
969  {
970  // Fall back on double precision when necessary since ExodusII
971  // doesn't seem to support long double
972  int comp_ws = std::min(sizeof(Real), sizeof(double));
973  int io_ws = std::min(sizeof(Real), sizeof(double));
974 
975  ex_id = exII::ex_create(filename.c_str(), EX_CLOBBER, &comp_ws, &io_ws);
976 
977  EX_CHECK_ERR(ex_id, "Error creating ExodusII mesh file.");
978 
979  if (verbose)
980  libMesh::out << "File created successfully." << std::endl;
981  }
982 
983  opened_for_writing = true;
984  current_filename = filename;
985 }
int libMesh::ExodusII_IO_Helper::get_block_id ( int  index)

Get the block number for the given block index.

Definition at line 437 of file exodusII_io_helper.C.

References block_ids.

Referenced by libMesh::ExodusII_IO::read(), and write_element_values().

438 {
439  libmesh_assert_less (static_cast<unsigned int>(index), block_ids.size());
440 
441  return block_ids[index];
442 }
std::string libMesh::ExodusII_IO_Helper::get_block_name ( int  index)

Get the block name for the given block index if supplied in the mesh file. Otherwise an empty string is returned.

Definition at line 446 of file exodusII_io_helper.C.

References block_ids, and id_to_block_names.

Referenced by libMesh::ExodusII_IO::read().

447 {
448  libmesh_assert_less (static_cast<unsigned int>(index), block_ids.size());
449 
450  return id_to_block_names[block_ids[index]];
451 }
const char * libMesh::ExodusII_IO_Helper::get_elem_type ( ) const
Returns
the current element type. Note: the default behavior is for this value to be in all capital letters, e.g. HEX27.

Definition at line 266 of file exodusII_io_helper.C.

References elem_type.

Referenced by libMesh::ExodusII_IO::read().

267 {
268  return &elem_type[0];
269 }
int libMesh::ExodusII_IO_Helper::get_node_set_id ( int  index)

Get the node set id for the given node set index.

Definition at line 473 of file exodusII_io_helper.C.

References nodeset_ids.

474 {
475  libmesh_assert_less (static_cast<unsigned int>(index), nodeset_ids.size());
476 
477  return nodeset_ids[index];
478 }
std::string libMesh::ExodusII_IO_Helper::get_node_set_name ( int  index)

Get the node set name for the given node set index if supplied in the mesh file. Otherwise an empty string is returned.

Definition at line 482 of file exodusII_io_helper.C.

References id_to_ns_names, and nodeset_ids.

Referenced by libMesh::ExodusII_IO::read().

483 {
484  libmesh_assert_less (static_cast<unsigned int>(index), nodeset_ids.size());
485 
486  return id_to_ns_names[nodeset_ids[index]];
487 }
int libMesh::ExodusII_IO_Helper::get_side_set_id ( int  index)

Get the side set id for the given side set index.

Definition at line 455 of file exodusII_io_helper.C.

References ss_ids.

Referenced by libMesh::ExodusII_IO::read().

456 {
457  libmesh_assert_less (static_cast<unsigned int>(index), ss_ids.size());
458 
459  return ss_ids[index];
460 }
std::string libMesh::ExodusII_IO_Helper::get_side_set_name ( int  index)

Get the side set name for the given side set index if supplied in the mesh file. Otherwise an empty string is returned.

Definition at line 464 of file exodusII_io_helper.C.

References id_to_ss_names, and ss_ids.

Referenced by libMesh::ExodusII_IO::read().

465 {
466  libmesh_assert_less (static_cast<unsigned int>(index), ss_ids.size());
467 
468  return id_to_ss_names[ss_ids[index]];
469 }
void libMesh::ExodusII_IO_Helper::initialize ( std::string  title,
const MeshBase mesh 
)
virtual

Initializes the Exodus file

Reimplemented in libMesh::Nemesis_IO_Helper.

Definition at line 1051 of file exodusII_io_helper.C.

References _run_only_on_proc0, _use_mesh_dimension_instead_of_spatial_dimension, libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::MeshBase::boundary_info, end, libMesh::err, ex_err, ex_id, ex_put_init(), libMesh::DofObject::id(), libMesh::MeshBase::mesh_dimension(), libMesh::MeshBase::n_active_elem(), libMesh::MeshBase::n_elem(), libMesh::MeshBase::n_nodes(), num_dim, num_elem, num_elem_blk, num_node_sets, num_nodes, num_side_sets, libMesh::ParallelObject::processor_id(), libMesh::MeshBase::spatial_dimension(), and libMesh::Elem::subdomain_id().

Referenced by libMesh::ExodusII_IO::write(), and libMesh::ExodusII_IO::write_nodal_data_common().

1052 {
1053  // n_active_elem() is a parallel_only function
1054  unsigned int n_active_elem = mesh.n_active_elem();
1055 
1056  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1057  return;
1058 
1060  num_dim = mesh.mesh_dimension();
1061  else
1062  num_dim = mesh.spatial_dimension();
1063 
1064  num_nodes = mesh.n_nodes();
1065  num_elem = mesh.n_elem();
1066 
1067  std::vector<boundary_id_type> unique_side_boundaries;
1068  std::vector<boundary_id_type> unique_node_boundaries;
1069 
1070  mesh.boundary_info->build_side_boundary_ids(unique_side_boundaries);
1071  mesh.boundary_info->build_node_boundary_ids(unique_node_boundaries);
1072 
1073  num_side_sets = unique_side_boundaries.size();
1074  num_node_sets = unique_node_boundaries.size();
1075 
1076  //loop through element and map between block and element vector
1077  std::map<subdomain_id_type, std::vector<unsigned int> > subdomain_map;
1078 
1079  MeshBase::const_element_iterator it = mesh.active_elements_begin();
1080  const MeshBase::const_element_iterator end = mesh.active_elements_end();
1081  for (; it!=end; ++it)
1082  {
1083  const Elem * elem = *it;
1084  subdomain_id_type cur_subdomain = elem->subdomain_id();
1085 
1086  subdomain_map[cur_subdomain].push_back(elem->id());
1087  }
1088  num_elem_blk = subdomain_map.size();
1089 
1090  if (str_title.size() > MAX_LINE_LENGTH)
1091  {
1092  libMesh::err << "Warning, Exodus files cannot have titles longer than "
1093  << MAX_LINE_LENGTH
1094  << " characters. Your title will be truncated."
1095  << std::endl;
1096  str_title.resize(MAX_LINE_LENGTH);
1097  }
1098 
1100  str_title.c_str(),
1101  num_dim,
1102  num_nodes,
1103  n_active_elem,
1104  num_elem_blk,
1105  num_node_sets,
1106  num_side_sets);
1107 
1108  EX_CHECK_ERR(ex_err, "Error initializing new Exodus file.");
1109 }
void libMesh::ExodusII_IO_Helper::initialize_discontinuous ( std::string  title,
const MeshBase mesh 
)

Initializes the Exodus file

Definition at line 990 of file exodusII_io_helper.C.

References _run_only_on_proc0, _use_mesh_dimension_instead_of_spatial_dimension, libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::MeshBase::boundary_info, end, libMesh::err, ex_err, ex_id, ex_put_init(), libMesh::DofObject::id(), libMesh::MeshBase::mesh_dimension(), libMesh::MeshBase::n_elem(), num_dim, num_elem, num_elem_blk, num_node_sets, num_nodes, num_side_sets, libMesh::ParallelObject::processor_id(), libMesh::MeshBase::spatial_dimension(), and libMesh::Elem::subdomain_id().

Referenced by libMesh::ExodusII_IO::write_nodal_data_common().

991 {
992  if ((_run_only_on_proc0) && (this->processor_id() != 0))
993  return;
994 
996  num_dim = mesh.mesh_dimension();
997  else
998  num_dim = mesh.spatial_dimension();
999 
1000  MeshBase::const_element_iterator it = mesh.active_elements_begin();
1001  const MeshBase::const_element_iterator end = mesh.active_elements_end();
1002  for (; it!=end; ++it)
1003  num_nodes += (*it)->n_nodes();
1004 
1005  num_elem = mesh.n_elem();
1006 
1007  std::vector<boundary_id_type> unique_side_boundaries;
1008  std::vector<boundary_id_type> unique_node_boundaries;
1009 
1010  mesh.boundary_info->build_side_boundary_ids(unique_side_boundaries);
1011  mesh.boundary_info->build_node_boundary_ids(unique_node_boundaries);
1012 
1013  num_side_sets = unique_side_boundaries.size();
1014  num_node_sets = unique_node_boundaries.size();
1015 
1016  //loop through element and map between block and element vector
1017  std::map<subdomain_id_type, std::vector<unsigned int> > subdomain_map;
1018 
1019  for (it=mesh.active_elements_begin(); it!=end; ++it)
1020  {
1021  const Elem * elem = *it;
1022  subdomain_id_type cur_subdomain = elem->subdomain_id();
1023 
1024  subdomain_map[cur_subdomain].push_back(elem->id());
1025  }
1026  num_elem_blk = subdomain_map.size();
1027 
1028  if (str_title.size() > MAX_LINE_LENGTH)
1029  {
1030  libMesh::err << "Warning, Exodus files cannot have titles longer than "
1031  << MAX_LINE_LENGTH
1032  << " characters. Your title will be truncated."
1033  << std::endl;
1034  str_title.resize(MAX_LINE_LENGTH);
1035  }
1036 
1038  str_title.c_str(),
1039  num_dim,
1040  num_nodes,
1041  num_elem,
1042  num_elem_blk,
1043  num_node_sets,
1044  num_side_sets);
1045 
1046  EX_CHECK_ERR(ex_err, "Error initializing new Exodus file.");
1047 }
void libMesh::ExodusII_IO_Helper::initialize_element_variables ( std::vector< std::string >  names)

Sets up the nodal variables

Definition at line 1553 of file exodusII_io_helper.C.

References _elem_vars_initialized, _run_only_on_proc0, check_existing_vars(), elem_var_names, ELEMENTAL, ex_err, ex_id, ex_put_elem_var_tab(), num_elem_blk, num_elem_vars, libMesh::ParallelObject::processor_id(), and write_var_names().

Referenced by libMesh::ExodusII_IO::write_element_data().

1554 {
1555  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1556  return;
1557 
1558  // Quick return if there are no element variables to write
1559  if (names.size() == 0)
1560  return;
1561 
1562  // Quick return if we have already called this function
1564  return;
1565 
1566  // Be sure that variables in the file match what we are asking for
1567  if (num_elem_vars > 0)
1568  {
1569  this->check_existing_vars(ELEMENTAL, names, this->elem_var_names);
1570  return;
1571  }
1572 
1573  // Set the flag so we can skip this stuff on subsequent calls to
1574  // initialize_element_variables()
1575  _elem_vars_initialized = true;
1576 
1577  this->write_var_names(ELEMENTAL, names);
1578 
1579  // Form the element variable truth table and send to Exodus.
1580  // This tells which variables are written to which blocks,
1581  // and can dramatically speed up writing element variables
1582  //
1583  // We really should initialize all entries in the truth table to 0
1584  // and then loop over all subdomains, setting their entries to 1
1585  // if a given variable exists on that subdomain. However,
1586  // we don't have that information, and the element variables
1587  // passed to us are padded with zeroes for the blocks where
1588  // they aren't defined. To be consistent with that, fill
1589  // the truth table with ones.
1590  std::vector<int> truth_tab(num_elem_blk*num_elem_vars, 1);
1592  num_elem_blk,
1593  num_elem_vars,
1594  &truth_tab[0]);
1595  EX_CHECK_ERR(ex_err, "Error writing element truth table.");
1596 }
void libMesh::ExodusII_IO_Helper::initialize_global_variables ( std::vector< std::string >  names)

Sets up the global variables

Definition at line 1628 of file exodusII_io_helper.C.

References _global_vars_initialized, _run_only_on_proc0, check_existing_vars(), GLOBAL, global_var_names, num_global_vars, libMesh::ParallelObject::processor_id(), and write_var_names().

Referenced by libMesh::Nemesis_IO::write_global_data(), and libMesh::ExodusII_IO::write_global_data().

1629 {
1630  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1631  return;
1632 
1633  // Quick return if there are no global variables to write
1634  if (names.size() == 0)
1635  return;
1636 
1638  return;
1639 
1640  // Be sure that variables in the file match what we are asking for
1641  if (num_global_vars > 0)
1642  {
1643  this->check_existing_vars(GLOBAL, names, this->global_var_names);
1644  return;
1645  }
1646 
1647  _global_vars_initialized = true;
1648 
1649  this->write_var_names(GLOBAL, names);
1650 }
void libMesh::ExodusII_IO_Helper::initialize_nodal_variables ( std::vector< std::string >  names)

Sets up the nodal variables

Definition at line 1600 of file exodusII_io_helper.C.

References _nodal_vars_initialized, _run_only_on_proc0, check_existing_vars(), NODAL, nodal_var_names, num_nodal_vars, libMesh::ParallelObject::processor_id(), and write_var_names().

Referenced by libMesh::Nemesis_IO::write_nodal_data(), and libMesh::ExodusII_IO::write_nodal_data_common().

1601 {
1602  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1603  return;
1604 
1605  // Quick return if there are no nodal variables to write
1606  if (names.size() == 0)
1607  return;
1608 
1609  // Quick return if we have already called this function
1611  return;
1612 
1613  // Be sure that variables in the file match what we are asking for
1614  if (num_nodal_vars > 0)
1615  {
1616  this->check_existing_vars(NODAL, names, this->nodal_var_names);
1617  return;
1618  }
1619 
1620  // Set the flag so we can skip the rest of this function on subsequent calls.
1621  _nodal_vars_initialized = true;
1622 
1623  this->write_var_names(NODAL, names);
1624 }
int libMesh::ExodusII_IO_Helper::inquire ( int  req_info,
std::string  error_msg = "" 
)

Generic inquiry, returns the value

Definition at line 703 of file exodusII_io_helper.C.

References ex_err, ex_id, and ex_inquire().

Referenced by read_num_time_steps(), read_sideset_info(), and write_information_records().

704 {
705  int ret_int = 0;
706  char ret_char = 0;
707  float ret_float = 0.;
708 
710  req_info_in,
711  &ret_int,
712  &ret_float,
713  &ret_char);
714 
715  EX_CHECK_ERR(ex_err, error_msg);
716 
717  return ret_int;
718 }
void libMesh::ExodusII_IO_Helper::message ( const std::string  msg)

Prints the message defined in msg. Can be turned off if verbosity is set to 0.

Definition at line 273 of file exodusII_io_helper.C.

References libMesh::out, and verbose.

Referenced by close(), read_block_info(), read_elem_in_block(), read_elem_num_map(), read_header(), read_node_num_map(), read_nodes(), read_nodeset(), read_nodeset_info(), read_sideset(), and read_sideset_info().

274 {
275  if (verbose) libMesh::out << msg << std::endl;
276 }
void libMesh::ExodusII_IO_Helper::message ( const std::string  msg,
int  i 
)

Prints the message defined in msg, and appends the number i to the end of the message. Useful for printing messages in loops. Can be turned off if verbosity is set to 0.

Definition at line 280 of file exodusII_io_helper.C.

References libMesh::out, and verbose.

281 {
282  if (verbose) libMesh::out << msg << i << "." << std::endl;
283 }
processor_id_type libMesh::ParallelObject::n_processors ( ) const
inlineinherited
Returns
the number of processors in the group.

Definition at line 92 of file parallel_object.h.

References libMesh::ParallelObject::_communicator, and libMesh::Parallel::Communicator::size().

Referenced by libMesh::ParmetisPartitioner::_do_repartition(), libMesh::ParallelMesh::add_elem(), libMesh::ParallelMesh::add_node(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::ParmetisPartitioner::assign_partitioning(), libMesh::ParallelMesh::assign_unique_ids(), libMesh::AztecLinearSolver< T >::AztecLinearSolver(), libMesh::ParallelMesh::clear(), libMesh::Nemesis_IO_Helper::compute_border_node_ids(), libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::UnstructuredMesh::create_pid_mesh(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::EnsightIO::EnsightIO(), libMesh::MeshBase::get_info(), libMesh::EquationSystems::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::ParmetisPartitioner::initialize(), libMesh::Nemesis_IO_Helper::initialize(), libMesh::MeshTools::libmesh_assert_valid_dof_ids(), libMesh::MeshTools::libmesh_assert_valid_procids< Elem >(), libMesh::MeshTools::libmesh_assert_valid_procids< Node >(), libMesh::MeshTools::libmesh_assert_valid_refinement_flags(), libMesh::MeshBase::n_active_elem_on_proc(), libMesh::MeshBase::n_elem_on_proc(), libMesh::MeshBase::n_nodes_on_proc(), libMesh::Partitioner::partition(), libMesh::MeshBase::partition(), libMesh::Partitioner::partition_unpartitioned_elements(), libMesh::PetscLinearSolver< T >::PetscLinearSolver(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::MeshTools::processor_bounding_box(), libMesh::System::project_vector(), libMesh::Nemesis_IO::read(), libMesh::CheckpointIO::read(), libMesh::UnstructuredMesh::read(), libMesh::System::read_parallel_data(), libMesh::System::read_SCALAR_dofs(), libMesh::System::read_serialized_blocked_dof_objects(), libMesh::System::read_serialized_vector(), libMesh::Partitioner::repartition(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::BoundaryInfo::sync(), libMesh::ParallelMesh::update_parallel_id_counts(), libMesh::CheckpointIO::write(), libMesh::GMVIO::write_binary(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::System::write_parallel_data(), libMesh::System::write_SCALAR_dofs(), libMesh::XdrIO::write_serialized_bcs(), libMesh::System::write_serialized_blocked_dof_objects(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), and libMesh::XdrIO::write_serialized_nodesets().

93  { return libmesh_cast_int<processor_id_type>(_communicator.size()); }
void libMesh::ExodusII_IO_Helper::open ( const char *  filename,
bool  read_only 
)

Opens an ExodusII mesh file named filename. If read_only==true, the file will be opened with the EX_READ flag, otherwise it will be opened with the EX_WRITE flag.

Definition at line 287 of file exodusII_io_helper.C.

References current_filename, ex_id, opened_for_reading, opened_for_writing, libMesh::out, libMesh::Real, and verbose.

Referenced by libMesh::ExodusII_IO::read(), libMesh::Nemesis_IO::read(), libMesh::Nemesis_IO::write_nodal_data(), and libMesh::ExodusII_IO::write_nodal_data_common().

288 {
289  // Version of Exodus you are using
290  float ex_version = 0.;
291 
292  // Word size in bytes of the floating point variables used in the
293  // application program (0, 4, or 8)
294  int comp_ws = sizeof(Real);
295 
296  // Word size in bytes of the floating point data as they are stored
297  // in the ExodusII file. "If this argument is 0, the word size of the
298  // floating point data already stored in the file is returned"
299  int io_ws = 0;
300 
301  ex_id = exII::ex_open(filename,
302  read_only ? EX_READ : EX_WRITE,
303  &comp_ws,
304  &io_ws,
305  &ex_version);
306 
307  std::string err_msg = std::string("Error opening ExodusII mesh file: ") + std::string(filename);
308  EX_CHECK_ERR(ex_id, err_msg);
309  if (verbose) libMesh::out << "File opened successfully." << std::endl;
310 
311  if (read_only)
312  opened_for_reading = true;
313  else
314  opened_for_writing = true;
315 
316  current_filename = std::string(filename);
317 }
void libMesh::ExodusII_IO_Helper::print_header ( )

Prints the ExodusII mesh file header, which includes the mesh title, the number of nodes, number of elements, mesh dimension, number of sidesets, and number of nodesets.

Definition at line 351 of file exodusII_io_helper.C.

References num_dim, num_elem, num_elem_blk, num_node_sets, num_nodes, num_side_sets, libMesh::out, title, and verbose.

Referenced by libMesh::ExodusII_IO::read(), and libMesh::Nemesis_IO::read().

352 {
353  if (verbose)
354  libMesh::out << "Title: \t" << &title[0] << std::endl
355  << "Mesh Dimension: \t" << num_dim << std::endl
356  << "Number of Nodes: \t" << num_nodes << std::endl
357  << "Number of elements: \t" << num_elem << std::endl
358  << "Number of elt blocks: \t" << num_elem_blk << std::endl
359  << "Number of node sets: \t" << num_node_sets << std::endl
360  << "Number of side sets: \t" << num_side_sets << std::endl;
361 }
void libMesh::ExodusII_IO_Helper::print_nodes ( std::ostream &  out = libMesh::out)

Prints the nodal information, by default to libMesh::out.

Definition at line 402 of file exodusII_io_helper.C.

References num_nodes, x, y, and z.

403 {
404  for (int i=0; i<num_nodes; i++)
405  out_stream << "(" << x[i] << ", " << y[i] << ", " << z[i] << ")" << std::endl;
406 }
processor_id_type libMesh::ParallelObject::processor_id ( ) const
inlineinherited
Returns
the rank of this processor in the group.

Definition at line 98 of file parallel_object.h.

References libMesh::ParallelObject::_communicator, and libMesh::Parallel::Communicator::rank().

Referenced by libMesh::MetisPartitioner::_do_partition(), libMesh::EquationSystems::_read_impl(), libMesh::SerialMesh::active_local_elements_begin(), libMesh::ParallelMesh::active_local_elements_begin(), libMesh::SerialMesh::active_local_elements_end(), libMesh::ParallelMesh::active_local_elements_end(), libMesh::SerialMesh::active_local_subdomain_elements_begin(), libMesh::ParallelMesh::active_local_subdomain_elements_begin(), libMesh::SerialMesh::active_local_subdomain_elements_end(), libMesh::ParallelMesh::active_local_subdomain_elements_end(), libMesh::SerialMesh::active_not_local_elements_begin(), libMesh::ParallelMesh::active_not_local_elements_begin(), libMesh::SerialMesh::active_not_local_elements_end(), libMesh::ParallelMesh::active_not_local_elements_end(), libMesh::ParallelMesh::add_elem(), libMesh::DofMap::add_neighbors_to_send_list(), libMesh::ParallelMesh::add_node(), libMesh::UnstructuredMesh::all_second_order(), libMesh::ParmetisPartitioner::assign_partitioning(), libMesh::EquationSystems::build_discontinuous_solution_vector(), libMesh::Nemesis_IO_Helper::build_element_and_node_maps(), libMesh::ParmetisPartitioner::build_graph(), libMesh::InfElemBuilder::build_inf_elem(), libMesh::DofMap::build_sparsity(), libMesh::ParallelMesh::clear(), close(), libMesh::Nemesis_IO_Helper::compute_border_node_ids(), libMesh::Nemesis_IO_Helper::compute_communication_map_parameters(), libMesh::Nemesis_IO_Helper::compute_internal_and_border_elems_and_internal_nodes(), libMesh::Nemesis_IO_Helper::compute_node_communication_maps(), libMesh::Nemesis_IO_Helper::compute_num_global_elem_blocks(), libMesh::Nemesis_IO_Helper::compute_num_global_nodesets(), libMesh::Nemesis_IO_Helper::compute_num_global_sidesets(), libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), create(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::DofMap::end_dof(), libMesh::DofMap::end_old_dof(), libMesh::EnsightIO::EnsightIO(), libMesh::UnstructuredMesh::find_neighbors(), libMesh::DofMap::first_dof(), libMesh::DofMap::first_old_dof(), libMesh::Nemesis_IO_Helper::get_cmap_params(), libMesh::Nemesis_IO_Helper::get_eb_info_global(), libMesh::Nemesis_IO_Helper::get_elem_cmap(), libMesh::Nemesis_IO_Helper::get_elem_map(), libMesh::MeshBase::get_info(), libMesh::Nemesis_IO_Helper::get_init_global(), libMesh::Nemesis_IO_Helper::get_init_info(), libMesh::Nemesis_IO_Helper::get_loadbal_param(), libMesh::Nemesis_IO_Helper::get_node_cmap(), libMesh::Nemesis_IO_Helper::get_node_map(), libMesh::Nemesis_IO_Helper::get_ns_param_global(), libMesh::Nemesis_IO_Helper::get_ss_param_global(), libMesh::MeshFunction::gradient(), libMesh::MeshFunction::hessian(), libMesh::SystemSubsetBySubdomain::init(), libMesh::ParmetisPartitioner::initialize(), initialize(), initialize_discontinuous(), initialize_element_variables(), initialize_global_variables(), initialize_nodal_variables(), libMesh::SparsityPattern::Build::join(), libMesh::DofMap::last_dof(), libMesh::MeshTools::libmesh_assert_valid_procids< Elem >(), libMesh::MeshTools::libmesh_assert_valid_procids< Node >(), libMesh::SerialMesh::local_elements_begin(), libMesh::ParallelMesh::local_elements_begin(), libMesh::SerialMesh::local_elements_end(), libMesh::ParallelMesh::local_elements_end(), libMesh::SerialMesh::local_level_elements_begin(), libMesh::ParallelMesh::local_level_elements_begin(), libMesh::SerialMesh::local_level_elements_end(), libMesh::ParallelMesh::local_level_elements_end(), libMesh::SerialMesh::local_nodes_begin(), libMesh::ParallelMesh::local_nodes_begin(), libMesh::SerialMesh::local_nodes_end(), libMesh::ParallelMesh::local_nodes_end(), libMesh::SerialMesh::local_not_level_elements_begin(), libMesh::ParallelMesh::local_not_level_elements_begin(), libMesh::SerialMesh::local_not_level_elements_end(), libMesh::ParallelMesh::local_not_level_elements_end(), libMesh::MeshRefinement::make_coarsening_compatible(), libMesh::MeshBase::n_active_local_elem(), libMesh::BoundaryInfo::n_boundary_conds(), libMesh::BoundaryInfo::n_edge_conds(), libMesh::DofMap::n_local_dofs(), libMesh::System::n_local_dofs(), libMesh::MeshBase::n_local_elem(), libMesh::MeshBase::n_local_nodes(), libMesh::BoundaryInfo::n_nodeset_conds(), libMesh::SerialMesh::not_local_elements_begin(), libMesh::ParallelMesh::not_local_elements_begin(), libMesh::SerialMesh::not_local_elements_end(), libMesh::ParallelMesh::not_local_elements_end(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::SparsityPattern::Build::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::MeshFunction::operator()(), libMesh::ParallelMesh::ParallelMesh(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::System::project_vector(), libMesh::Nemesis_IO_Helper::put_cmap_params(), libMesh::Nemesis_IO_Helper::put_elem_cmap(), libMesh::Nemesis_IO_Helper::put_elem_map(), libMesh::Nemesis_IO_Helper::put_loadbal_param(), libMesh::Nemesis_IO_Helper::put_node_cmap(), libMesh::Nemesis_IO_Helper::put_node_map(), libMesh::Nemesis_IO::read(), libMesh::CheckpointIO::read(), libMesh::XdrIO::read(), libMesh::UnstructuredMesh::read(), libMesh::CheckpointIO::read_connectivity(), read_elem_num_map(), libMesh::System::read_header(), libMesh::System::read_legacy_data(), read_node_num_map(), libMesh::System::read_parallel_data(), libMesh::System::read_SCALAR_dofs(), libMesh::XdrIO::read_serialized_bc_names(), libMesh::XdrIO::read_serialized_bcs(), libMesh::System::read_serialized_blocked_dof_objects(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::System::read_serialized_data(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::System::read_serialized_vector(), libMesh::System::read_serialized_vectors(), libMesh::MeshData::read_xdr(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::BoundaryInfo::sync(), libMesh::MeshTools::total_weight(), libMesh::ParallelMesh::update_parallel_id_counts(), libMesh::MeshTools::weight(), libMesh::ExodusII_IO::write(), libMesh::CheckpointIO::write(), libMesh::XdrIO::write(), libMesh::UnstructuredMesh::write(), libMesh::EquationSystems::write(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::ExodusII_IO::write_element_data(), write_element_values(), write_elements(), write_elements_discontinuous(), libMesh::ExodusII_IO::write_global_data(), write_global_values(), libMesh::System::write_header(), libMesh::ExodusII_IO::write_information_records(), write_information_records(), write_nodal_coordinates(), write_nodal_coordinates_discontinuous(), libMesh::UCDIO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data_discontinuous(), write_nodal_values(), write_nodesets(), libMesh::Nemesis_IO_Helper::write_nodesets(), libMesh::System::write_parallel_data(), libMesh::System::write_SCALAR_dofs(), libMesh::XdrIO::write_serialized_bc_names(), libMesh::XdrIO::write_serialized_bcs(), libMesh::System::write_serialized_blocked_dof_objects(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::System::write_serialized_data(), libMesh::XdrIO::write_serialized_nodes(), libMesh::XdrIO::write_serialized_nodesets(), libMesh::XdrIO::write_serialized_subdomain_names(), libMesh::System::write_serialized_vector(), libMesh::System::write_serialized_vectors(), write_sidesets(), libMesh::Nemesis_IO_Helper::write_sidesets(), libMesh::ExodusII_IO::write_timestep(), and write_timestep().

99  { return libmesh_cast_int<processor_id_type>(_communicator.rank()); }
void libMesh::ExodusII_IO_Helper::read_block_info ( )

Reads information for all of the blocks in the ExodusII mesh file.

Definition at line 410 of file exodusII_io_helper.C.

References block_ids, EX_ELEM_BLOCK, ex_err, ex_get_elem_blk_ids(), ex_get_name(), ex_id, id_to_block_names, message(), and num_elem_blk.

Referenced by libMesh::ExodusII_IO::read(), libMesh::Nemesis_IO::read(), and libMesh::ExodusII_IO::write_nodal_data_common().

411 {
412  block_ids.resize(num_elem_blk);
413  // Get all element block IDs.
415  block_ids.empty() ? NULL : &block_ids[0]);
416  // Usually, there is only one
417  // block since there is only
418  // one type of element.
419  // However, there could be more.
420 
421  EX_CHECK_ERR(ex_err, "Error getting block IDs.");
422  message("All block IDs retrieved successfully.");
423 
424  char name_buffer[MAX_STR_LENGTH+1];
425  for (int i=0; i<num_elem_blk; ++i)
426  {
428  block_ids[i], name_buffer);
429  EX_CHECK_ERR(ex_err, "Error getting block name.");
430  id_to_block_names[block_ids[i]] = name_buffer;
431  }
432  message("All block names retrieved successfully.");
433 }
void libMesh::ExodusII_IO_Helper::read_elem_in_block ( int  block)

Reads all of the element connectivity for block block in the ExodusII mesh file.

Definition at line 492 of file exodusII_io_helper.C.

References block_ids, connect, elem_type, ex_err, ex_get_elem_block(), ex_get_elem_conn(), ex_id, message(), num_attr, num_elem_this_blk, num_nodes_per_elem, libMesh::out, and verbose.

Referenced by libMesh::ExodusII_IO::read(), and libMesh::Nemesis_IO::read().

493 {
494  libmesh_assert_less (static_cast<unsigned int>(block), block_ids.size());
495 
497  block_ids[block],
498  &elem_type[0],
501  &num_attr);
502  if (verbose)
503  libMesh::out << "Reading a block of " << num_elem_this_blk
504  << " " << &elem_type[0] << "(s)"
505  << " having " << num_nodes_per_elem
506  << " nodes per element." << std::endl;
507 
508  EX_CHECK_ERR(ex_err, "Error getting block info.");
509  message("Info retrieved successfully for block: ", block);
510 
511 
512 
513  // Read in the connectivity of the elements of this block,
514  // watching out for the case where we actually have no
515  // elements in this block (possible with parallel files)
517 
518  if (!connect.empty())
519  {
520  ex_err = exII::ex_get_elem_conn(ex_id,
521  block_ids[block],
522  &connect[0]);
523 
524  EX_CHECK_ERR(ex_err, "Error reading block connectivity.");
525  message("Connectivity retrieved successfully for block: ", block);
526  }
527 }
void libMesh::ExodusII_IO_Helper::read_elem_num_map ( )

Reads the optional node_num_map from the ExodusII mesh file.

Definition at line 532 of file exodusII_io_helper.C.

References elem_num_map, ex_err, ex_get_elem_num_map(), ex_id, message(), std::min(), num_elem, libMesh::out, libMesh::ParallelObject::processor_id(), and verbose.

Referenced by libMesh::ExodusII_IO::read(), and libMesh::Nemesis_IO::read().

533 {
534  elem_num_map.resize(num_elem);
535 
537  elem_num_map.empty() ? NULL : &elem_num_map[0]);
538 
539  EX_CHECK_ERR(ex_err, "Error retrieving element number map.");
540  message("Element numbering map retrieved successfully.");
541 
542 
543  if (verbose)
544  {
545  libMesh::out << "[" << this->processor_id() << "] elem_num_map[i] = ";
546  for (unsigned int i=0; i<static_cast<unsigned int>(std::min(10, num_elem-1)); ++i)
547  libMesh::out << elem_num_map[i] << ", ";
548  libMesh::out << "... " << elem_num_map.back() << std::endl;
549  }
550 }
void libMesh::ExodusII_IO_Helper::read_elemental_var_values ( std::string  elemental_var_name,
int  time_step 
)

Reads elemental values for the variable 'elemental_var_name' at the specified timestep into the 'elem_var_values' array.

Definition at line 900 of file exodusII_io_helper.C.

References block_ids, elem_var_names, elem_var_values, ELEMENTAL, libMesh::err, ex_err, ex_get_elem_block(), ex_get_elem_var(), ex_id, num_elem, num_elem_blk, and read_var_names().

Referenced by libMesh::ExodusII_IO::copy_elemental_solution().

901 {
902  // CAUTION: this assumes that libMesh element numbering is identical to exodus block-by-block element numbering
903  // There is no way how to get the whole elemental field from the exodus file, so we have to go block by block
904 
905  elem_var_values.resize(num_elem);
906 
907  this->read_var_names(ELEMENTAL);
908 
909  // See if we can find the variable we are looking for
910  unsigned int var_index = 0;
911  bool found = false;
912 
913  // Do a linear search for nodal_var_name in nodal_var_names
914  for (; var_index<elem_var_names.size(); ++var_index)
915  {
916  found = (elem_var_names[var_index] == elemental_var_name);
917  if (found)
918  break;
919  }
920 
921  if (!found)
922  {
923  libMesh::err << "Unable to locate variable named: " << elemental_var_name << std::endl;
924  libMesh::err << "Available variables: " << std::endl;
925  for (unsigned int i=0; i<elem_var_names.size(); ++i)
926  libMesh::err << elem_var_names[i] << std::endl;
927 
928  libmesh_error();
929  }
930 
931  unsigned int ex_el_num = 0;
932  for (unsigned int i=0; i<static_cast<unsigned int>(num_elem_blk); i++)
933  {
934  int n_blk_elems = 0;
936  block_ids[i],
937  NULL,
938  &n_blk_elems,
939  NULL,
940  NULL);
941  EX_CHECK_ERR(ex_err, "Error getting number of elements in block.");
942 
943  std::vector<Real> block_elem_var_values(num_elem);
945  time_step,
946  var_index+1,
947  block_ids[i],
948  n_blk_elems,
949  &block_elem_var_values[0]);
950  EX_CHECK_ERR(ex_err, "Error getting elemental values.");
951 
952  for (unsigned int j=0; j<static_cast<unsigned int>(n_blk_elems); j++)
953  {
954  elem_var_values[ex_el_num] = block_elem_var_values[j];
955  ex_el_num++;
956  }
957  }
958 }
void libMesh::ExodusII_IO_Helper::read_header ( )

Reads an ExodusII mesh file header.

Definition at line 321 of file exodusII_io_helper.C.

References ex_err, ex_get_init(), ex_get_var_param(), ex_id, message(), num_dim, num_elem, num_elem_blk, num_elem_vars, num_global_vars, num_nodal_vars, num_node_sets, num_nodes, num_side_sets, read_num_time_steps(), and title.

Referenced by libMesh::ExodusII_IO::read(), libMesh::Nemesis_IO::read(), libMesh::Nemesis_IO::write_nodal_data(), and libMesh::ExodusII_IO::write_nodal_data_common().

322 {
324  &title[0],
325  &num_dim,
326  &num_nodes,
327  &num_elem,
328  &num_elem_blk,
329  &num_node_sets,
330  &num_side_sets);
331 
332  EX_CHECK_ERR(ex_err, "Error retrieving header info.");
333 
334  this->read_num_time_steps();
335 
337  EX_CHECK_ERR(ex_err, "Error reading number of nodal variables.");
338 
340  EX_CHECK_ERR(ex_err, "Error reading number of elemental variables.");
341 
343  EX_CHECK_ERR(ex_err, "Error reading number of global variables.");
344 
345  message("Exodus header info retrieved successfully.");
346 }
void libMesh::ExodusII_IO_Helper::read_nodal_var_values ( std::string  nodal_var_name,
int  time_step 
)

Reads the nodal values for the variable 'nodal_var_name' at the specified time into the 'nodal_var_values' array.

Definition at line 745 of file exodusII_io_helper.C.

References libMesh::err, ex_err, ex_get_nodal_var(), ex_id, NODAL, nodal_var_names, nodal_var_values, num_nodes, and read_var_names().

Referenced by libMesh::ExodusII_IO::copy_nodal_solution().

746 {
747  // Read the nodal variable names from file, so we can see if we have the one we're looking for
748  this->read_var_names(NODAL);
749 
750  // See if we can find the variable we are looking for
751  unsigned int var_index = 0;
752  bool found = false;
753 
754  // Do a linear search for nodal_var_name in nodal_var_names
755  for (; var_index<nodal_var_names.size(); ++var_index)
756  {
757  found = (nodal_var_names[var_index] == nodal_var_name);
758  if (found)
759  break;
760  }
761 
762  if (!found)
763  {
764  libMesh::err << "Unable to locate variable named: " << nodal_var_name << std::endl;
765  libMesh::err << "Available variables: " << std::endl;
766  for (unsigned int i=0; i<nodal_var_names.size(); ++i)
767  libMesh::err << nodal_var_names[i] << std::endl;
768 
769  libmesh_error();
770  }
771 
772  // Allocate enough space to store the nodal variable values
773  nodal_var_values.resize(num_nodes);
774 
775  // Call the Exodus API to read the nodal variable values
777  time_step,
778  var_index+1,
779  num_nodes,
780  &nodal_var_values[0]);
781  EX_CHECK_ERR(ex_err, "Error reading nodal variable values!");
782 }
void libMesh::ExodusII_IO_Helper::read_node_num_map ( )

Reads the optional node_num_map from the ExodusII mesh file.

Definition at line 382 of file exodusII_io_helper.C.

References ex_err, ex_get_node_num_map(), ex_id, message(), std::min(), node_num_map, num_nodes, libMesh::out, libMesh::ParallelObject::processor_id(), and verbose.

Referenced by libMesh::ExodusII_IO::read(), and libMesh::Nemesis_IO::read().

383 {
384  node_num_map.resize(num_nodes);
385 
387  node_num_map.empty() ? NULL : &node_num_map[0]);
388 
389  EX_CHECK_ERR(ex_err, "Error retrieving nodal number map.");
390  message("Nodal numbering map retrieved successfully.");
391 
392  if (verbose)
393  {
394  libMesh::out << "[" << this->processor_id() << "] node_num_map[i] = ";
395  for (unsigned int i=0; i<static_cast<unsigned int>(std::min(10, num_nodes-1)); ++i)
396  libMesh::out << node_num_map[i] << ", ";
397  libMesh::out << "... " << node_num_map.back() << std::endl;
398  }
399 }
void libMesh::ExodusII_IO_Helper::read_nodes ( )

Reads the nodal data (x,y,z coordinates) from the ExodusII mesh file.

Definition at line 365 of file exodusII_io_helper.C.

References ex_err, ex_get_coord(), ex_id, message(), num_nodes, x, y, and z.

Referenced by libMesh::ExodusII_IO::read(), and libMesh::Nemesis_IO::read().

366 {
367  x.resize(num_nodes);
368  y.resize(num_nodes);
369  z.resize(num_nodes);
370 
372  static_cast<void*>(&x[0]),
373  static_cast<void*>(&y[0]),
374  static_cast<void*>(&z[0]));
375 
376  EX_CHECK_ERR(ex_err, "Error retrieving nodal data.");
377  message("Nodal data retrieved successfully.");
378 }
void libMesh::ExodusII_IO_Helper::read_nodeset ( int  id)

Reads information about nodeset id and inserts it into the global nodeset array at the position offset.

Definition at line 658 of file exodusII_io_helper.C.

References ex_err, ex_get_node_set(), ex_get_node_set_param(), ex_id, message(), node_list, nodeset_ids, num_node_df_per_set, and num_nodes_per_set.

Referenced by libMesh::ExodusII_IO::read(), and libMesh::Nemesis_IO::read().

659 {
660  libmesh_assert_less (static_cast<unsigned int>(id), nodeset_ids.size());
661  libmesh_assert_less (static_cast<unsigned int>(id), num_nodes_per_set.size());
662  libmesh_assert_less (static_cast<unsigned int>(id), num_node_df_per_set.size());
663 
665  nodeset_ids[id],
666  &num_nodes_per_set[id],
667  &num_node_df_per_set[id]);
668  EX_CHECK_ERR(ex_err, "Error retrieving nodeset parameters.");
669  message("Parameters retrieved successfully for nodeset: ", id);
670 
671  node_list.resize(num_nodes_per_set[id]);
672 
673  // Don't call ex_get_node_set unless there are actually nodes there to get.
674  // Exodus prints an annoying warning message in DEBUG mode otherwise...
675  if (num_nodes_per_set[id] > 0)
676  {
677  ex_err = exII::ex_get_node_set(ex_id,
678  nodeset_ids[id],
679  &node_list[0]);
680 
681  EX_CHECK_ERR(ex_err, "Error retrieving nodeset data.");
682  message("Data retrieved successfully for nodeset: ", id);
683  }
684 }
void libMesh::ExodusII_IO_Helper::read_nodeset_info ( )

Reads information about all of the nodesets in the ExodusII mesh file.

Definition at line 587 of file exodusII_io_helper.C.

References ex_err, ex_get_name(), ex_get_node_set_ids(), ex_id, EX_NODE_SET, id_to_ns_names, message(), nodeset_ids, num_node_df_per_set, num_node_sets, and num_nodes_per_set.

Referenced by libMesh::ExodusII_IO::read(), and libMesh::Nemesis_IO::read().

588 {
589  nodeset_ids.resize(num_node_sets);
590  if (num_node_sets > 0)
591  {
593  &nodeset_ids[0]);
594  EX_CHECK_ERR(ex_err, "Error retrieving nodeset information.");
595  message("All nodeset information retrieved successfully.");
596 
597  // Resize appropriate data structures -- only do this once outnode the loop
600  }
601 
602  char name_buffer[MAX_STR_LENGTH+1];
603  for (int i=0; i<num_node_sets; ++i)
604  {
606  nodeset_ids[i], name_buffer);
607  EX_CHECK_ERR(ex_err, "Error getting node set name.");
608  id_to_ns_names[nodeset_ids[i]] = name_buffer;
609  }
610  message("All node set names retrieved successfully.");
611 }
void libMesh::ExodusII_IO_Helper::read_num_time_steps ( )

Reads the number of timesteps currently stored in the Exodus file and stores it in the num_time_steps variable.

Definition at line 737 of file exodusII_io_helper.C.

References EX_INQ_TIME, inquire(), and num_time_steps.

Referenced by libMesh::ExodusII_IO::get_num_time_steps(), read_header(), and read_time_steps().

738 {
740  this->inquire(exII::EX_INQ_TIME, "Error retrieving number of time steps");
741 }
void libMesh::ExodusII_IO_Helper::read_sideset ( int  id,
int  offset 
)

Reads information about sideset id and inserts it into the global sideset array at the position offset.

Definition at line 615 of file exodusII_io_helper.C.

References elem_list, ex_err, ex_get_side_set(), ex_get_side_set_param(), ex_id, id_list, message(), num_df_per_set, num_sides_per_set, side_list, and ss_ids.

Referenced by libMesh::ExodusII_IO::read(), and libMesh::Nemesis_IO::read().

616 {
617  libmesh_assert_less (static_cast<unsigned int>(id), ss_ids.size());
618  libmesh_assert_less (static_cast<unsigned int>(id), num_sides_per_set.size());
619  libmesh_assert_less (static_cast<unsigned int>(id), num_df_per_set.size());
620  libmesh_assert_less_equal (static_cast<unsigned int>(offset), elem_list.size());
621  libmesh_assert_less_equal (static_cast<unsigned int>(offset), side_list.size());
622 
624  ss_ids[id],
625  &num_sides_per_set[id],
626  &num_df_per_set[id]);
627  EX_CHECK_ERR(ex_err, "Error retrieving sideset parameters.");
628  message("Parameters retrieved successfully for sideset: ", id);
629 
630 
631  // It's OK for offset==elem_list.size() as long as num_sides_per_set[id]==0
632  // because in that case we don't actually read anything...
633 #ifdef DEBUG
634  if (static_cast<unsigned int>(offset) == elem_list.size() ||
635  static_cast<unsigned int>(offset) == side_list.size() )
636  libmesh_assert_equal_to (num_sides_per_set[id], 0);
637 #endif
638 
639 
640  // Don't call ex_get_side_set unless there are actually sides there to get.
641  // Exodus prints an annoying warning in DEBUG mode otherwise...
642  if (num_sides_per_set[id] > 0)
643  {
644  ex_err = exII::ex_get_side_set(ex_id,
645  ss_ids[id],
646  &elem_list[offset],
647  &side_list[offset]);
648  EX_CHECK_ERR(ex_err, "Error retrieving sideset data.");
649  message("Data retrieved successfully for sideset: ", id);
650 
651  for (int i=0; i<num_sides_per_set[id]; i++)
652  id_list[i+offset] = ss_ids[id];
653  }
654 }
void libMesh::ExodusII_IO_Helper::read_sideset_info ( )

Reads information about all of the sidesets in the ExodusII mesh file.

Definition at line 554 of file exodusII_io_helper.C.

References elem_list, ex_err, ex_get_name(), ex_get_side_set_ids(), ex_id, EX_INQ_SS_ELEM_LEN, EX_SIDE_SET, id_list, id_to_ss_names, inquire(), message(), num_df_per_set, num_elem_all_sidesets, num_side_sets, num_sides_per_set, side_list, and ss_ids.

Referenced by libMesh::ExodusII_IO::read(), and libMesh::Nemesis_IO::read().

555 {
556  ss_ids.resize(num_side_sets);
557  if (num_side_sets > 0)
558  {
560  &ss_ids[0]);
561  EX_CHECK_ERR(ex_err, "Error retrieving sideset information.");
562  message("All sideset information retrieved successfully.");
563 
564  // Resize appropriate data structures -- only do this once outside the loop
567 
568  // Inquire about the length of the concatenated side sets element list
569  num_elem_all_sidesets = inquire(exII::EX_INQ_SS_ELEM_LEN, "Error retrieving length of the concatenated side sets element list!");
570 
574  }
575 
576  char name_buffer[MAX_STR_LENGTH+1];
577  for (int i=0; i<num_side_sets; ++i)
578  {
580  ss_ids[i], name_buffer);
581  EX_CHECK_ERR(ex_err, "Error getting side set name.");
582  id_to_ss_names[ss_ids[i]] = name_buffer;
583  }
584  message("All side set names retrieved successfully.");
585 }
void libMesh::ExodusII_IO_Helper::read_time_steps ( )

Reads and stores the timesteps in the 'time_steps' array.

Definition at line 722 of file exodusII_io_helper.C.

References ex_err, ex_get_all_times(), ex_id, num_time_steps, read_num_time_steps(), and time_steps.

Referenced by libMesh::ExodusII_IO::get_time_steps().

723 {
724  // Make sure we have an up-to-date count of the number of time steps in the file.
725  this->read_num_time_steps();
726 
727  if (num_time_steps > 0)
728  {
729  time_steps.resize(num_time_steps);
731  EX_CHECK_ERR(ex_err, "Error reading timesteps!");
732  }
733 }
void libMesh::ExodusII_IO_Helper::read_var_names ( ExodusVarType  type)
private

Definition at line 786 of file exodusII_io_helper.C.

References elem_var_names, ELEMENTAL, libMesh::err, GLOBAL, global_var_names, NODAL, nodal_var_names, num_elem_vars, num_global_vars, num_nodal_vars, and read_var_names_impl().

Referenced by check_existing_vars(), read_elemental_var_values(), and read_nodal_var_values().

787 {
788  switch (type)
789  {
790  case NODAL:
792  break;
793  case ELEMENTAL:
795  break;
796  case GLOBAL:
798  break;
799  default:
800  libMesh::err << "Unrecognized ExodusVarType " << type << std::endl;
801  libmesh_error();
802  }
803 }
void libMesh::ExodusII_IO_Helper::read_var_names_impl ( const char *  var_type,
int &  count,
std::vector< std::string > &  result 
)
private

read_var_names() dispatches to this function.

Definition at line 807 of file exodusII_io_helper.C.

References ex_err, ex_get_var_names(), ex_get_var_param(), ex_id, libMesh::ExodusII_IO_Helper::NamesData::get_char_star(), libMesh::ExodusII_IO_Helper::NamesData::get_char_star_star(), libMesh::out, and verbose.

Referenced by read_var_names().

810 {
811  // First read and store the number of names we have
812  ex_err = exII::ex_get_var_param(ex_id, var_type, &count);
813  EX_CHECK_ERR(ex_err, "Error reading number of variables.");
814 
815  // Second read the actual names and convert them into a format we can use
816  NamesData names_table(count, MAX_STR_LENGTH);
817 
819  var_type,
820  count,
821  names_table.get_char_star_star()
822  );
823  EX_CHECK_ERR(ex_err, "Error reading variable names!");
824 
825  if (verbose)
826  {
827  libMesh::out << "Read the variable(s) from the file:" << std::endl;
828  for (int i=0; i<count; i++)
829  libMesh::out << names_table.get_char_star(i) << std::endl;
830  }
831 
832  // Allocate enough space for our variable name strings.
833  result.resize(count);
834 
835  // Copy the char buffers into strings.
836  for (int i=0; i<count; i++)
837  result[i] = names_table.get_char_star(i); // calls string::op=(const char*)
838 }
void libMesh::ExodusII_IO_Helper::set_coordinate_offset ( Point  p)

Allows you to set a vector that is added to the coordinates of all of the nodes. Effectively, this "moves" the mesh to a particular position

Definition at line 1822 of file exodusII_io_helper.C.

References _coordinate_offset.

Referenced by libMesh::ExodusII_IO::set_coordinate_offset().

1823 {
1824  _coordinate_offset = p;
1825 }
void libMesh::ExodusII_IO_Helper::use_mesh_dimension_instead_of_spatial_dimension ( bool  val)

Sets the underlying value of the boolean flag _use_mesh_dimension_instead_of_spatial_dimension. By default, the value of this flag is false.

See the ExodusII_IO class documentation for a detailed description of this flag.

Definition at line 1817 of file exodusII_io_helper.C.

References _use_mesh_dimension_instead_of_spatial_dimension.

Referenced by libMesh::ExodusII_IO::use_mesh_dimension_instead_of_spatial_dimension().

1818 {
1820 }
void libMesh::ExodusII_IO_Helper::write_element_values ( const MeshBase mesh,
const std::vector< Number > &  values,
int  timestep 
)

Writes the vector of values to the element variables.

Definition at line 1697 of file exodusII_io_helper.C.

References _run_only_on_proc0, libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), data, end, ex_err, ex_id, ex_put_elem_var(), ex_update(), get_block_id(), libMesh::DofObject::id(), num_elem, libMesh::ParallelObject::processor_id(), and libMesh::Elem::subdomain_id().

Referenced by libMesh::ExodusII_IO::write_element_data().

1698 {
1699  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1700  return;
1701 
1702  // Loop over the element blocks and write the data one block at a time
1703  std::map<unsigned int, std::vector<unsigned int> > subdomain_map;
1704 
1705  const unsigned int num_vars = values.size() / num_elem;
1706 
1707  MeshBase::const_element_iterator mesh_it = mesh.active_elements_begin();
1708  const MeshBase::const_element_iterator end = mesh.active_elements_end();
1709 
1710  // loop through element and map between block and element vector
1711  for (; mesh_it!=end; ++mesh_it)
1712  {
1713  const Elem * elem = *mesh_it;
1714  subdomain_map[elem->subdomain_id()].push_back(elem->id());
1715  }
1716 
1717  // For each variable, create a 'data' array which holds all the elemental variable
1718  // values *for a given block* on this processor, then write that data vector to file
1719  // before moving onto the next block.
1720  for (unsigned int i=0; i<num_vars; ++i)
1721  {
1722  // The size of the subdomain map is the number of blocks.
1723  std::map<unsigned int, std::vector<unsigned int> >::iterator it = subdomain_map.begin();
1724 
1725  for (unsigned int j=0; it!=subdomain_map.end(); ++it, ++j)
1726  {
1727  const std::vector<unsigned int> & elem_nums = (*it).second;
1728  const unsigned int num_elems_this_block = elem_nums.size();
1729  std::vector<Number> data(num_elems_this_block);
1730 
1731  for (unsigned int k=0; k<num_elems_this_block; ++k)
1732  data[k] = values[i*num_elem + elem_nums[k]];
1733 
1735  timestep,
1736  i+1,
1737  this->get_block_id(j),
1738  num_elems_this_block,
1739  &data[0]);
1740  EX_CHECK_ERR(ex_err, "Error writing element values.");
1741  }
1742  }
1743 
1745  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
1746 }
void libMesh::ExodusII_IO_Helper::write_elements ( const MeshBase mesh)
virtual

Writes the elements contained in "mesh". FIXME: This only works for Meshes having a single type of element!

Reimplemented in libMesh::Nemesis_IO_Helper.

Definition at line 1210 of file exodusII_io_helper.C.

References _run_only_on_proc0, libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::ExodusII_IO_Helper::ElementMaps::assign_conversion(), block_ids, connect, libMesh::MeshBase::elem(), elem_num_map, end, libMesh::Utility::enum_to_string(), libMesh::err, EX_ELEM_BLOCK, ex_err, ex_id, ex_put_elem_block(), ex_put_elem_conn(), ex_put_elem_num_map(), ex_put_names(), libMesh::ExodusII_IO_Helper::Conversion::exodus_elem_type(), libMesh::ExodusII_IO_Helper::Conversion::get_canonical_type(), libMesh::ExodusII_IO_Helper::NamesData::get_char_star_star(), libMesh::ExodusII_IO_Helper::Conversion::get_inverse_node_map(), libMesh::DofObject::id(), libmesh_elem_num_to_exodus, libMesh::MeshBase::n_active_elem(), libMesh::Elem::n_nodes(), libMesh::Elem::node(), num_elem_blk, num_nodes_per_elem, libMesh::out, libMesh::ParallelObject::processor_id(), libMesh::ExodusII_IO_Helper::NamesData::push_back_entry(), libMesh::Elem::subdomain_id(), libMesh::MeshBase::subdomain_name(), libMesh::Elem::type(), and verbose.

Referenced by libMesh::ExodusII_IO::write(), and libMesh::ExodusII_IO::write_nodal_data_common().

1211 {
1212  // n_active_elem() is a parallel_only function
1213  unsigned int n_active_elem = mesh.n_active_elem();
1214 
1215  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1216  return;
1217 
1218  // Map from block ID to a vector of element IDs in that block. Element
1219  // IDs are now of type dof_id_type, subdomain IDs are of type subdomain_id_type.
1220  typedef std::map<subdomain_id_type, std::vector<dof_id_type> > subdomain_map_type;
1221  subdomain_map_type subdomain_map;
1222 
1223  MeshBase::const_element_iterator mesh_it = mesh.active_elements_begin();
1224  const MeshBase::const_element_iterator end = mesh.active_elements_end();
1225  //loop through element and map between block and element vector
1226  for (; mesh_it!=end; ++mesh_it)
1227  {
1228  const Elem * elem = *mesh_it;
1229  subdomain_map[ elem->subdomain_id() ].push_back(elem->id());
1230  }
1231 
1232  // element map vector
1233  num_elem_blk = subdomain_map.size();
1234  block_ids.resize(num_elem_blk);
1235  elem_num_map.resize(n_active_elem);
1236  std::vector<int>::iterator curr_elem_map_end = elem_num_map.begin();
1237 
1238  // Note: It appears that there is a bug in exodusII::ex_put_name where
1239  // the index returned from the ex_id_lkup is erronously used. For now
1240  // the work around is to use the alternative function ex_put_names, but
1241  // this function requires a char** datastructure.
1242  NamesData names_table(num_elem_blk, MAX_STR_LENGTH);
1243 
1244  // This counter is used to fill up the libmesh_elem_num_to_exodus map in the loop below.
1245  unsigned libmesh_elem_num_to_exodus_counter = 0;
1246 
1247  // counter indexes into the block_ids vector
1248  unsigned int counter = 0;
1249 
1250  for (subdomain_map_type::iterator it=subdomain_map.begin(); it!=subdomain_map.end(); ++it)
1251  {
1252  block_ids[counter] = (*it).first;
1253  names_table.push_back_entry(mesh.subdomain_name((*it).first));
1254 
1255  // Get a reference to a vector of element IDs for this subdomain.
1256  subdomain_map_type::mapped_type& tmp_vec = (*it).second;
1257 
1258  ExodusII_IO_Helper::ElementMaps em;
1259 
1260  //Use the first element in this block to get representative information.
1261  //Note that Exodus assumes all elements in a block are of the same type!
1262  //We are using that same assumption here!
1263  const ExodusII_IO_Helper::Conversion conv = em.assign_conversion(mesh.elem(tmp_vec[0])->type());
1264  num_nodes_per_elem = mesh.elem(tmp_vec[0])->n_nodes();
1265 
1267  (*it).first,
1268  conv.exodus_elem_type().c_str(),
1269  tmp_vec.size(),
1271  /*num_attr=*/0);
1272 
1273  EX_CHECK_ERR(ex_err, "Error writing element block.");
1274 
1275  connect.resize(tmp_vec.size()*num_nodes_per_elem);
1276 
1277  for (unsigned int i=0; i<tmp_vec.size(); i++)
1278  {
1279  unsigned int elem_id = tmp_vec[i];
1280  libmesh_elem_num_to_exodus[elem_id] = ++libmesh_elem_num_to_exodus_counter; // 1-based indexing for Exodus
1281 
1282  const Elem* elem = mesh.elem(elem_id);
1283 
1284  // We *might* be able to get away with writing mixed element
1285  // types which happen to have the same number of nodes, but
1286  // do we actually *want* to get away with that?
1287  // .) No visualization software would be able to handle it.
1288  // .) There'd be no way for us to read it back in reliably.
1289  // .) Even elements with the same number of nodes may have different connectivities (?)
1290 
1291  // This needs to be more than an assert so we don't fail
1292  // with a mysterious segfault while trying to write mixed
1293  // element meshes in optimized mode.
1294  if (elem->type() != conv.get_canonical_type())
1295  {
1296  libMesh::err << "Error: Exodus requires all elements with a given subdomain ID to be the same type.\n"
1297  << "Can't write both "
1298  << Utility::enum_to_string(elem->type())
1299  << " and "
1300  << Utility::enum_to_string(conv.get_canonical_type())
1301  << " in the same block!"
1302  << std::endl;
1303  libmesh_error();
1304  }
1305 
1306  for (unsigned int j=0; j<static_cast<unsigned int>(num_nodes_per_elem); ++j)
1307  {
1308  unsigned connect_index = (i*num_nodes_per_elem)+j;
1309  unsigned elem_node_index = conv.get_inverse_node_map(j); // inverse node map is for writing.
1310  if (verbose)
1311  {
1312  libMesh::out << "Exodus node index: " << j
1313  << "=LibMesh node index " << elem_node_index << std::endl;
1314  }
1315 
1316  // FIXME: We are hard-coding the 1-based node numbering assumption here.
1317  connect[connect_index] = elem->node(elem_node_index)+1;
1318  }
1319  }
1320  ex_err = exII::ex_put_elem_conn(ex_id, (*it).first, &connect[0]);
1321  EX_CHECK_ERR(ex_err, "Error writing element connectivities");
1322 
1323  // This transform command stores its result in a range that begins at the third argument,
1324  // so this command is adding values to the elem_num_map vector starting from curr_elem_map_end.
1325  curr_elem_map_end = std::transform(tmp_vec.begin(),
1326  tmp_vec.end(),
1327  curr_elem_map_end,
1328  std::bind2nd(std::plus<subdomain_map_type::mapped_type::value_type>(), 1)); // Adds one to each id to make a 1-based exodus file!
1329 
1330  // But if we don't want to add one, we just want to put the values
1331  // of tmp_vec into elem_map in the right location, we can use
1332  // std::copy().
1333  // curr_elem_map_end = std::copy(tmp_vec.begin(), tmp_vec.end(), curr_elem_map_end);
1334 
1335  counter++;
1336  }
1337 
1338  // write out the element number map that we created
1340  EX_CHECK_ERR(ex_err, "Error writing element map");
1341 
1342  // Write out the block names
1343  if (num_elem_blk > 0)
1344  {
1345  ex_err = exII::ex_put_names(ex_id, exII::EX_ELEM_BLOCK, names_table.get_char_star_star());
1346  EX_CHECK_ERR(ex_err, "Error writing element names");
1347  }
1348 }
void libMesh::ExodusII_IO_Helper::write_elements_discontinuous ( const MeshBase mesh)

Writes the elements contained in "mesh". FIXME: This only works for Meshes having a single type of element!

Definition at line 1353 of file exodusII_io_helper.C.

References _run_only_on_proc0, libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::ExodusII_IO_Helper::ElementMaps::assign_conversion(), connect, libMesh::MeshBase::elem(), elem_num_map, end, ex_err, ex_id, ex_put_elem_block(), ex_put_elem_conn(), ex_put_elem_num_map(), ex_update(), libMesh::ExodusII_IO_Helper::Conversion::exodus_elem_type(), libMesh::ExodusII_IO_Helper::Conversion::get_inverse_node_map(), libMesh::DofObject::id(), libmesh_elem_num_to_exodus, libMesh::Elem::n_nodes(), num_nodes_per_elem, libMesh::out, libMesh::ParallelObject::processor_id(), libMesh::Elem::subdomain_id(), libMesh::Elem::type(), and verbose.

Referenced by libMesh::ExodusII_IO::write_nodal_data_common().

1354 {
1355  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1356  return;
1357 
1358  typedef std::map<subdomain_id_type, std::vector<dof_id_type> > subdomain_map_type;
1359  subdomain_map_type subdomain_map;
1360 
1361  MeshBase::const_element_iterator mesh_it = mesh.active_elements_begin();
1362  const MeshBase::const_element_iterator end = mesh.active_elements_end();
1363  // loop through element and map between block and element vector
1364  for (; mesh_it!=end; ++mesh_it)
1365  {
1366  const Elem * elem = *mesh_it;
1367  subdomain_map[ elem->subdomain_id() ].push_back(elem->id());
1368  }
1369 
1370  // This counter is used to fill up the libmesh_elem_num_to_exodus map in the loop below.
1371  unsigned libmesh_elem_num_to_exodus_counter = 0;
1372 
1373  for (subdomain_map_type::iterator it=subdomain_map.begin(); it!=subdomain_map.end(); ++it)
1374  {
1375  subdomain_map_type::mapped_type& tmp_vec = (*it).second;
1376 
1377  ExodusII_IO_Helper::ElementMaps em;
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!
1382  const ExodusII_IO_Helper::Conversion conv = em.assign_conversion(mesh.elem(tmp_vec[0])->type());
1383  num_nodes_per_elem = mesh.elem(tmp_vec[0])->n_nodes();
1384 
1386  (*it).first,
1387  conv.exodus_elem_type().c_str(),
1388  tmp_vec.size(),
1390  /*num_attr=*/0);
1391 
1392  EX_CHECK_ERR(ex_err, "Error writing element block.");
1393 
1394  connect.resize(tmp_vec.size()*num_nodes_per_elem);
1395 
1396  for (unsigned int i=0; i<tmp_vec.size(); i++)
1397  {
1398  unsigned int elem_id = tmp_vec[i];
1399  libmesh_elem_num_to_exodus[elem_id] = ++libmesh_elem_num_to_exodus_counter; // 1-based indexing for Exodus
1400 
1401  for (unsigned int j=0; j<static_cast<unsigned int>(num_nodes_per_elem); j++)
1402  {
1403  const unsigned int connect_index = (i*num_nodes_per_elem)+j;
1404  const unsigned int elem_node_index = conv.get_inverse_node_map(j); // Inverse node map is for writing
1405  if (verbose)
1406  {
1407  libMesh::out << "Exodus node index: " << j
1408  << "=LibMesh node index " << elem_node_index << std::endl;
1409  }
1410  connect[connect_index] = i*num_nodes_per_elem+elem_node_index+1;
1411  }
1412  }
1413  ex_err = exII::ex_put_elem_conn(ex_id, (*it).first, &connect[0]);
1414  EX_CHECK_ERR(ex_err, "Error writing element connectivities");
1415 
1416  // Create space in elem_num_map
1417  elem_num_map.resize(tmp_vec.size());
1418 
1419  // copy the contents of tmp_vec into the elem_map vector
1420  std::copy(tmp_vec.begin(), tmp_vec.end(), elem_num_map.begin());
1421 
1422  // And write to file
1424  EX_CHECK_ERR(ex_err, "Error writing element map");
1425  }
1426 
1428  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
1429 }
void libMesh::ExodusII_IO_Helper::write_global_values ( const std::vector< Number > &  values,
int  timestep 
)

Writes the vector of global variables.

Definition at line 1803 of file exodusII_io_helper.C.

References _run_only_on_proc0, ex_err, ex_id, ex_put_glob_vars(), ex_update(), num_global_vars, and libMesh::ParallelObject::processor_id().

Referenced by libMesh::Nemesis_IO::write_global_data(), and libMesh::ExodusII_IO::write_global_data().

1804 {
1805  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1806  return;
1807 
1808  ex_err = exII::ex_put_glob_vars(ex_id, timestep, num_global_vars, &values[0]);
1809  EX_CHECK_ERR(ex_err, "Error writing global values.");
1810 
1812  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
1813 }
void libMesh::ExodusII_IO_Helper::write_information_records ( const std::vector< std::string > &  records)

Writes the vector of information records.

Definition at line 1764 of file exodusII_io_helper.C.

References _run_only_on_proc0, libMesh::err, ex_err, ex_id, EX_INQ_INFO, ex_put_info(), ex_update(), libMesh::ExodusII_IO_Helper::NamesData::get_char_star_star(), inquire(), libMesh::ParallelObject::processor_id(), and libMesh::ExodusII_IO_Helper::NamesData::push_back_entry().

Referenced by libMesh::Nemesis_IO::write_information_records(), and libMesh::ExodusII_IO::write_information_records().

1765 {
1766  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1767  return;
1768 
1769  // There may already be information records in the file (for
1770  // example, if we're appending) and in that case, according to the
1771  // Exodus documentation, writing more information records is not
1772  // supported.
1773  int num_info = inquire(exII::EX_INQ_INFO, "Error retrieving the number of information records from file!");
1774  if (num_info > 0)
1775  {
1776  libMesh::err << "Warning! The Exodus file already contains information records.\n"
1777  << "Exodus does not support writing additional records in this situation."
1778  << std::endl;
1779  return;
1780  }
1781 
1782  int num_records = records.size();
1783 
1784  if (num_records > 0)
1785  {
1786  NamesData info(num_records, MAX_LINE_LENGTH);
1787 
1788  // If an entry is longer than MAX_LINE_LENGTH characters it's not an error, we just
1789  // write the first MAX_LINE_LENGTH characters to the file.
1790  for (unsigned i=0; i<records.size(); ++i)
1791  info.push_back_entry(records[i]);
1792 
1793  ex_err = exII::ex_put_info(ex_id, num_records, info.get_char_star_star());
1794  EX_CHECK_ERR(ex_err, "Error writing global values.");
1795 
1797  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
1798  }
1799 }
void libMesh::ExodusII_IO_Helper::write_nodal_coordinates ( const MeshBase mesh)
virtual

Writes the nodal coordinates contained in "mesh"

Reimplemented in libMesh::Nemesis_IO_Helper.

Definition at line 1113 of file exodusII_io_helper.C.

References _coordinate_offset, _run_only_on_proc0, end, ex_err, ex_id, ex_put_coord(), ex_put_node_num_map(), libMesh::DofObject::id(), node_num_map, libMesh::MeshBase::nodes_begin(), libMesh::MeshBase::nodes_end(), num_nodes, libMesh::ParallelObject::processor_id(), x, y, and z.

Referenced by libMesh::ExodusII_IO::write(), and libMesh::ExodusII_IO::write_nodal_data_common().

1114 {
1115  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1116  return;
1117 
1118  // Make room in the coordinate vectors
1119  x.resize(num_nodes);
1120  y.resize(num_nodes);
1121  z.resize(num_nodes);
1122 
1123  // And in the node_num_map - since the nodes aren't organized in
1124  // blocks, libmesh will always write out the identity map
1125  // here... unless there has been some refinement and coarsening. If
1126  // the nodes aren't renumbered after refinement and coarsening,
1127  // there may be 'holes' in the numbering, so we write out the node
1128  // map just to be on the safe side.
1129  node_num_map.resize(num_nodes);
1130 
1131  {
1132  MeshBase::const_node_iterator it = mesh.nodes_begin();
1133  const MeshBase::const_node_iterator end = mesh.nodes_end();
1134  for (unsigned i = 0; it != end; ++it, ++i)
1135  {
1136  const Node* node = *it;
1137 
1138  x[i] = (*node)(0) + _coordinate_offset(0);
1139 
1140 #if LIBMESH_DIM > 1
1141  y[i]=(*node)(1) + _coordinate_offset(1);
1142 #else
1143  y[i]=0.;
1144 #endif
1145 #if LIBMESH_DIM > 2
1146  z[i]=(*node)(2) + _coordinate_offset(2);
1147 #else
1148  z[i]=0.;
1149 #endif
1150 
1151  // Fill in node_num_map entry with the proper (1-based) node id
1152  node_num_map[i] = node->id() + 1;
1153  }
1154  }
1155 
1157  x.empty() ? NULL : &x[0],
1158  y.empty() ? NULL : &y[0],
1159  z.empty() ? NULL : &z[0]);
1160 
1161  EX_CHECK_ERR(ex_err, "Error writing coordinates to Exodus file.");
1162 
1163  // Also write the (1-based) node_num_map to the file.
1165  EX_CHECK_ERR(ex_err, "Error writing node_num_map");
1166 }
void libMesh::ExodusII_IO_Helper::write_nodal_coordinates_discontinuous ( const MeshBase mesh)

Writes the nodal coordinates contained in "mesh"

Definition at line 1170 of file exodusII_io_helper.C.

References _run_only_on_proc0, libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), end, ex_err, ex_id, ex_put_coord(), num_nodes, libMesh::ParallelObject::processor_id(), x, y, and z.

Referenced by libMesh::ExodusII_IO::write_nodal_data_common().

1171 {
1172  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1173  return;
1174 
1175  x.resize(num_nodes);
1176  y.resize(num_nodes);
1177  z.resize(num_nodes);
1178 
1179  MeshBase::const_element_iterator it = mesh.active_elements_begin();
1180  const MeshBase::const_element_iterator end = mesh.active_elements_end();
1181 
1182  unsigned int i = 0;
1183  for (; it!=end; ++it)
1184  for (unsigned int n=0; n<(*it)->n_nodes(); n++)
1185  {
1186  x[i]=(*it)->point(n)(0);
1187 #if LIBMESH_DIM > 1
1188  y[i]=(*it)->point(n)(1);
1189 #else
1190  y[i]=0.;
1191 #endif
1192 #if LIBMESH_DIM > 2
1193  z[i]=(*it)->point(n)(2);
1194 #else
1195  z[i]=0.;
1196 #endif
1197  i++;
1198  }
1199 
1201  x.empty() ? NULL : &x[0],
1202  y.empty() ? NULL : &y[0],
1203  z.empty() ? NULL : &z[0]);
1204 
1205  EX_CHECK_ERR(ex_err, "Error writing coordinates to Exodus file.");
1206 }
void libMesh::ExodusII_IO_Helper::write_nodal_values ( int  var_id,
const std::vector< Number > &  values,
int  timestep 
)

Writes the vector of values to a nodal variable.

Definition at line 1750 of file exodusII_io_helper.C.

References _run_only_on_proc0, ex_err, ex_id, ex_put_nodal_var(), ex_update(), num_nodes, and libMesh::ParallelObject::processor_id().

Referenced by libMesh::ExodusII_IO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data_discontinuous(), and libMesh::Nemesis_IO_Helper::write_nodal_solution().

1751 {
1752  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1753  return;
1754 
1755  ex_err = exII::ex_put_nodal_var(ex_id, timestep, var_id, num_nodes, &values[0]);
1756  EX_CHECK_ERR(ex_err, "Error writing nodal values.");
1757 
1759  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
1760 }
void libMesh::ExodusII_IO_Helper::write_nodesets ( const MeshBase mesh)
virtual

Writes the nodesets contained in "mesh"

Reimplemented in libMesh::Nemesis_IO_Helper.

Definition at line 1505 of file exodusII_io_helper.C.

References _run_only_on_proc0, libMesh::MeshBase::boundary_info, ex_err, ex_id, EX_NODE_SET, ex_put_names(), ex_put_node_set(), ex_put_node_set_param(), libMesh::ParallelObject::processor_id(), and libMesh::ExodusII_IO_Helper::NamesData::push_back_entry().

Referenced by libMesh::ExodusII_IO::write(), and libMesh::ExodusII_IO::write_nodal_data_common().

1506 {
1507  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1508  return;
1509 
1510  std::vector< dof_id_type > nl;
1511  std::vector< boundary_id_type > il;
1512 
1513  mesh.boundary_info->build_node_list(nl, il);
1514 
1515  // Maps from nodeset id to the nodes
1516  std::map<boundary_id_type, std::vector<int> > node;
1517 
1518  // Accumulate the vectors to pass into ex_put_node_set
1519  for (unsigned int i=0; i<nl.size(); i++)
1520  node[il[i]].push_back(nl[i]+1);
1521 
1522  std::vector<boundary_id_type> node_boundary_ids;
1523  mesh.boundary_info->build_node_boundary_ids(node_boundary_ids);
1524 
1525  // Write out the nodeset names, but only if there is something to write
1526  if (node_boundary_ids.size() > 0)
1527  {
1528  NamesData names_table(node_boundary_ids.size(), MAX_STR_LENGTH);
1529 
1530  for (unsigned int i=0; i<node_boundary_ids.size(); i++)
1531  {
1532  int nodeset_id = node_boundary_ids[i];
1533 
1534  int actual_id = nodeset_id;
1535 
1536  names_table.push_back_entry(mesh.boundary_info->nodeset_name(nodeset_id));
1537 
1538  ex_err = exII::ex_put_node_set_param(ex_id, actual_id, node[nodeset_id].size(), 0);
1539  EX_CHECK_ERR(ex_err, "Error writing nodeset parameters");
1540 
1541  ex_err = exII::ex_put_node_set(ex_id, actual_id, &node[nodeset_id][0]);
1542  EX_CHECK_ERR(ex_err, "Error writing nodesets");
1543  }
1544 
1545  // Write out the nodeset names
1546  ex_err = exII::ex_put_names(ex_id, exII::EX_NODE_SET, names_table.get_char_star_star());
1547  EX_CHECK_ERR(ex_err, "Error writing nodeset names");
1548  }
1549 }
void libMesh::ExodusII_IO_Helper::write_sidesets ( const MeshBase mesh)
virtual

Writes the sidesets contained in "mesh"

We need to build up active elements if AMR is enabled and add them to the exodus sidesets instead of the potentially inactive "parent" elements

Reimplemented in libMesh::Nemesis_IO_Helper.

Definition at line 1433 of file exodusII_io_helper.C.

References _run_only_on_proc0, libMesh::Elem::active_family_tree_by_side(), libMesh::ExodusII_IO_Helper::ElementMaps::assign_conversion(), libMesh::MeshBase::boundary_info, libMesh::MeshBase::elem(), ex_err, ex_id, ex_put_names(), ex_put_side_set(), ex_put_side_set_param(), EX_SIDE_SET, libMesh::ExodusII_IO_Helper::Conversion::get_inverse_side_map(), libmesh_elem_num_to_exodus, libMesh::ParallelObject::processor_id(), libMesh::ExodusII_IO_Helper::NamesData::push_back_entry(), and side.

Referenced by libMesh::ExodusII_IO::write(), and libMesh::ExodusII_IO::write_nodal_data_common().

1434 {
1435  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1436  return;
1437 
1438  ExodusII_IO_Helper::ElementMaps em;
1439 
1440  std::vector< dof_id_type > el;
1441  std::vector< unsigned short int > sl;
1442  std::vector< boundary_id_type > il;
1443 
1444  mesh.boundary_info->build_side_list(el, sl, il);
1445 
1446  // Maps from sideset id to the element and sides
1447  std::map<int, std::vector<int> > elem;
1448  std::map<int, std::vector<int> > side;
1449 
1450  // Accumulate the vectors to pass into ex_put_side_set
1451  for (unsigned int i=0; i<el.size(); i++)
1452  {
1453  std::vector<const Elem *> family;
1454 #ifdef LIBMESH_ENABLE_AMR
1455 
1459  mesh.elem(el[i])->active_family_tree_by_side(family, sl[i], false);
1460 #else
1461  family.push_back(mesh.elem(el[i]));
1462 #endif
1463 
1464  for (unsigned int j=0; j<family.size(); ++j)
1465  {
1466  const ExodusII_IO_Helper::Conversion conv = em.assign_conversion(mesh.elem(family[j]->id())->type());
1467 
1468  // Use the libmesh to exodus datastructure map to get the proper sideset IDs
1469  // The datastructure contains the "collapsed" contiguous ids
1470  elem[il[i]].push_back(libmesh_elem_num_to_exodus[family[j]->id()]);
1471  side[il[i]].push_back(conv.get_inverse_side_map(sl[i]));
1472  }
1473  }
1474 
1475  std::vector<boundary_id_type> side_boundary_ids;
1476  mesh.boundary_info->build_side_boundary_ids(side_boundary_ids);
1477 
1478  // Write out the sideset names, but only if there is something to write
1479  if (side_boundary_ids.size() > 0)
1480  {
1481  NamesData names_table(side_boundary_ids.size(), MAX_STR_LENGTH);
1482 
1483  for (unsigned int i=0; i<side_boundary_ids.size(); i++)
1484  {
1485  int ss_id = side_boundary_ids[i];
1486 
1487  int actual_id = ss_id;
1488 
1489  names_table.push_back_entry(mesh.boundary_info->sideset_name(ss_id));
1490 
1491  ex_err = exII::ex_put_side_set_param(ex_id, actual_id, elem[ss_id].size(), 0);
1492  EX_CHECK_ERR(ex_err, "Error writing sideset parameters");
1493 
1494  ex_err = exII::ex_put_side_set(ex_id, actual_id, &elem[ss_id][0], &side[ss_id][0]);
1495  EX_CHECK_ERR(ex_err, "Error writing sidesets");
1496  }
1497 
1498  ex_err = exII::ex_put_names(ex_id, exII::EX_SIDE_SET, names_table.get_char_star_star());
1499  EX_CHECK_ERR(ex_err, "Error writing sideset names");
1500  }
1501 }
void libMesh::ExodusII_IO_Helper::write_timestep ( int  timestep,
Real  time 
)

Writes the time for the timestep

Definition at line 1683 of file exodusII_io_helper.C.

References _run_only_on_proc0, ex_err, ex_id, ex_put_time(), ex_update(), and libMesh::ParallelObject::processor_id().

Referenced by libMesh::Nemesis_IO::write_timestep(), and libMesh::ExodusII_IO::write_timestep().

1684 {
1685  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1686  return;
1687 
1688  ex_err = exII::ex_put_time(ex_id, timestep, &time);
1689  EX_CHECK_ERR(ex_err, "Error writing timestep.");
1690 
1692  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
1693 }
void libMesh::ExodusII_IO_Helper::write_var_names ( ExodusVarType  type,
std::vector< std::string > &  names 
)
private

Wraps calls to exII::ex_put_var_names() and exII::ex_put_var_param(). The enumeration controls whether nodal, elemental, or global variable names are read and which class members are filled in.

Definition at line 843 of file exodusII_io_helper.C.

References ELEMENTAL, libMesh::err, GLOBAL, NODAL, num_elem_vars, num_global_vars, num_nodal_vars, and write_var_names_impl().

Referenced by initialize_element_variables(), initialize_global_variables(), and initialize_nodal_variables().

844 {
845  switch (type)
846  {
847  case NODAL:
848  this->write_var_names_impl("n", num_nodal_vars, names);
849  break;
850  case ELEMENTAL:
851  this->write_var_names_impl("e", num_elem_vars, names);
852  break;
853  case GLOBAL:
854  this->write_var_names_impl("g", num_global_vars, names);
855  break;
856  default:
857  libMesh::err << "Unrecognized ExodusVarType " << type << std::endl;
858  libmesh_error();
859  }
860 }
void libMesh::ExodusII_IO_Helper::write_var_names_impl ( const char *  var_type,
int &  count,
std::vector< std::string > &  names 
)
private

write_var_names() dispatches to this function.

Definition at line 864 of file exodusII_io_helper.C.

References ex_err, ex_id, ex_put_var_names(), ex_put_var_param(), libMesh::out, libMesh::ExodusII_IO_Helper::NamesData::push_back_entry(), and verbose.

Referenced by write_var_names().

865 {
866  // Update the count variable so that it's available to other parts of the class.
867  count = names.size();
868 
869  // Write that number of variables to the file.
870  ex_err = exII::ex_put_var_param(ex_id, var_type, count);
871  EX_CHECK_ERR(ex_err, "Error setting number of vars.");
872 
873  if (names.size() > 0)
874  {
875  NamesData names_table(names.size(), MAX_STR_LENGTH);
876 
877  // Store the input names in the format required by Exodus.
878  for (unsigned i=0; i<names.size(); ++i)
879  names_table.push_back_entry(names[i]);
880 
881  if (verbose)
882  {
883  libMesh::out << "Writing variable name(s) to file: " << std::endl;
884  for (unsigned i=0; i<names.size(); ++i)
885  libMesh::out << names_table.get_char_star(i) << std::endl;
886  }
887 
889  var_type,
890  names.size(),
891  names_table.get_char_star_star()
892  );
893 
894  EX_CHECK_ERR(ex_err, "Error writing variable names.");
895  }
896 }

Member Data Documentation

Point libMesh::ExodusII_IO_Helper::_coordinate_offset
protected

Definition at line 546 of file exodusII_io_helper.h.

Referenced by set_coordinate_offset(), and write_nodal_coordinates().

bool libMesh::ExodusII_IO_Helper::_elem_vars_initialized
protected

Definition at line 532 of file exodusII_io_helper.h.

Referenced by initialize_element_variables().

bool libMesh::ExodusII_IO_Helper::_global_vars_initialized
protected

Definition at line 535 of file exodusII_io_helper.h.

Referenced by initialize_global_variables().

bool libMesh::ExodusII_IO_Helper::_nodal_vars_initialized
protected

Definition at line 538 of file exodusII_io_helper.h.

Referenced by initialize_nodal_variables().

bool libMesh::ExodusII_IO_Helper::_use_mesh_dimension_instead_of_spatial_dimension
protected
std::vector<int> libMesh::ExodusII_IO_Helper::connect
std::string libMesh::ExodusII_IO_Helper::current_filename
std::vector<int> libMesh::ExodusII_IO_Helper::elem_list
std::vector<int> libMesh::ExodusII_IO_Helper::elem_num_map
std::vector<char> libMesh::ExodusII_IO_Helper::elem_type
std::vector<std::string> libMesh::ExodusII_IO_Helper::elem_var_names
std::vector<Real> libMesh::ExodusII_IO_Helper::elem_var_values
int libMesh::ExodusII_IO_Helper::ex_id

Definition at line 369 of file exodusII_io_helper.h.

Referenced by close(), create(), libMesh::Nemesis_IO_Helper::create(), libMesh::Nemesis_IO_Helper::get_cmap_params(), libMesh::Nemesis_IO_Helper::get_eb_info_global(), libMesh::Nemesis_IO_Helper::get_elem_cmap(), libMesh::Nemesis_IO_Helper::get_elem_map(), libMesh::Nemesis_IO_Helper::get_init_global(), libMesh::Nemesis_IO_Helper::get_init_info(), libMesh::Nemesis_IO_Helper::get_loadbal_param(), libMesh::Nemesis_IO_Helper::get_node_cmap(), libMesh::Nemesis_IO_Helper::get_node_map(), libMesh::Nemesis_IO_Helper::get_ns_param_global(), libMesh::Nemesis_IO_Helper::get_ss_param_global(), initialize(), initialize_discontinuous(), initialize_element_variables(), inquire(), open(), libMesh::Nemesis_IO_Helper::put_cmap_params(), libMesh::Nemesis_IO_Helper::put_eb_info_global(), libMesh::Nemesis_IO_Helper::put_elem_cmap(), libMesh::Nemesis_IO_Helper::put_elem_map(), libMesh::Nemesis_IO_Helper::put_init_global(), libMesh::Nemesis_IO_Helper::put_init_info(), libMesh::Nemesis_IO_Helper::put_loadbal_param(), libMesh::Nemesis_IO_Helper::put_n_coord(), libMesh::Nemesis_IO_Helper::put_node_cmap(), libMesh::Nemesis_IO_Helper::put_node_map(), libMesh::Nemesis_IO_Helper::put_ns_param_global(), libMesh::Nemesis_IO_Helper::put_ss_param_global(), read_block_info(), read_elem_in_block(), read_elem_num_map(), read_elemental_var_values(), read_header(), read_nodal_var_values(), read_node_num_map(), read_nodes(), read_nodeset(), read_nodeset_info(), read_sideset(), read_sideset_info(), read_time_steps(), read_var_names_impl(), libMesh::Nemesis_IO::write(), write_element_values(), write_elements(), libMesh::Nemesis_IO_Helper::write_elements(), write_elements_discontinuous(), libMesh::Nemesis_IO_Helper::write_exodus_initialization_info(), write_global_values(), write_information_records(), write_nodal_coordinates(), libMesh::Nemesis_IO_Helper::write_nodal_coordinates(), write_nodal_coordinates_discontinuous(), write_nodal_values(), write_nodesets(), libMesh::Nemesis_IO_Helper::write_nodesets(), write_sidesets(), libMesh::Nemesis_IO_Helper::write_sidesets(), write_timestep(), write_var_names_impl(), and libMesh::Nemesis_IO_Helper::~Nemesis_IO_Helper().

std::vector<int> libMesh::ExodusII_IO_Helper::exodus_elem_num_to_libmesh
std::vector<int> libMesh::ExodusII_IO_Helper::exodus_node_num_to_libmesh
std::vector<std::string> libMesh::ExodusII_IO_Helper::global_var_names

Definition at line 499 of file exodusII_io_helper.h.

Referenced by initialize_global_variables(), and read_var_names().

std::vector<int> libMesh::ExodusII_IO_Helper::id_list
std::map<int, std::string> libMesh::ExodusII_IO_Helper::id_to_block_names

Definition at line 502 of file exodusII_io_helper.h.

Referenced by get_block_name(), and read_block_info().

std::map<int, std::string> libMesh::ExodusII_IO_Helper::id_to_ns_names

Definition at line 504 of file exodusII_io_helper.h.

Referenced by get_node_set_name(), and read_nodeset_info().

std::map<int, std::string> libMesh::ExodusII_IO_Helper::id_to_ss_names

Definition at line 503 of file exodusII_io_helper.h.

Referenced by get_side_set_name(), and read_sideset_info().

std::vector<std::string> libMesh::ExodusII_IO_Helper::nodal_var_names
std::vector<Real> libMesh::ExodusII_IO_Helper::nodal_var_values
std::vector<int> libMesh::ExodusII_IO_Helper::node_list
std::vector<int> libMesh::ExodusII_IO_Helper::node_num_map
std::vector<int> libMesh::ExodusII_IO_Helper::nodeset_ids
int libMesh::ExodusII_IO_Helper::num_attr

Definition at line 402 of file exodusII_io_helper.h.

Referenced by read_elem_in_block().

std::vector<int> libMesh::ExodusII_IO_Helper::num_df_per_set

Definition at line 426 of file exodusII_io_helper.h.

Referenced by read_sideset(), and read_sideset_info().

int libMesh::ExodusII_IO_Helper::num_dim
int libMesh::ExodusII_IO_Helper::num_elem_all_sidesets

Definition at line 405 of file exodusII_io_helper.h.

Referenced by libMesh::Nemesis_IO::read(), and read_sideset_info().

int libMesh::ExodusII_IO_Helper::num_elem_this_blk
int libMesh::ExodusII_IO_Helper::num_elem_vars
int libMesh::ExodusII_IO_Helper::num_global_vars
int libMesh::ExodusII_IO_Helper::num_nodal_vars
std::vector<int> libMesh::ExodusII_IO_Helper::num_node_df_per_set

Definition at line 429 of file exodusII_io_helper.h.

Referenced by read_nodeset(), and read_nodeset_info().

std::vector<int> libMesh::ExodusII_IO_Helper::num_nodes_per_set

Definition at line 423 of file exodusII_io_helper.h.

Referenced by read_nodeset(), and read_nodeset_info().

std::vector<int> libMesh::ExodusII_IO_Helper::num_sides_per_set
int libMesh::ExodusII_IO_Helper::num_time_steps
std::vector<int> libMesh::ExodusII_IO_Helper::side_list
std::vector<int> libMesh::ExodusII_IO_Helper::ss_ids
std::vector<Real> libMesh::ExodusII_IO_Helper::time_steps

Definition at line 478 of file exodusII_io_helper.h.

Referenced by libMesh::ExodusII_IO::get_time_steps(), and read_time_steps().

std::vector<char> libMesh::ExodusII_IO_Helper::title

Definition at line 459 of file exodusII_io_helper.h.

Referenced by ExodusII_IO_Helper(), print_header(), and read_header().

bool libMesh::ExodusII_IO_Helper::verbose

Definition at line 507 of file exodusII_io_helper.h.

Referenced by libMesh::Nemesis_IO_Helper::build_element_and_node_maps(), libMesh::Nemesis_IO_Helper::compute_border_node_ids(), libMesh::Nemesis_IO_Helper::compute_communication_map_parameters(), libMesh::Nemesis_IO_Helper::compute_internal_and_border_elems_and_internal_nodes(), libMesh::Nemesis_IO_Helper::compute_node_communication_maps(), libMesh::Nemesis_IO_Helper::compute_num_global_elem_blocks(), libMesh::Nemesis_IO_Helper::compute_num_global_nodesets(), libMesh::Nemesis_IO_Helper::compute_num_global_sidesets(), libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), create(), libMesh::Nemesis_IO_Helper::create(), libMesh::Nemesis_IO_Helper::get_cmap_params(), libMesh::Nemesis_IO_Helper::get_eb_info_global(), libMesh::Nemesis_IO_Helper::get_elem_cmap(), libMesh::Nemesis_IO_Helper::get_elem_map(), libMesh::Nemesis_IO_Helper::get_init_global(), libMesh::Nemesis_IO_Helper::get_init_info(), libMesh::Nemesis_IO_Helper::get_loadbal_param(), libMesh::Nemesis_IO_Helper::get_node_cmap(), libMesh::Nemesis_IO_Helper::get_node_map(), libMesh::Nemesis_IO_Helper::get_ns_param_global(), libMesh::Nemesis_IO_Helper::get_ss_param_global(), message(), open(), print_header(), libMesh::Nemesis_IO_Helper::put_node_cmap(), read_elem_in_block(), read_elem_num_map(), read_node_num_map(), read_var_names_impl(), libMesh::ExodusII_IO::verbose(), libMesh::Nemesis_IO::verbose(), write_elements(), write_elements_discontinuous(), libMesh::Nemesis_IO_Helper::write_nodesets(), libMesh::Nemesis_IO_Helper::write_sidesets(), and write_var_names_impl().


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

Site Created By: libMesh Developers
Last modified: February 07 2014 16:58:00 UTC

Hosted By:
SourceForge.net Logo