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 Types

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

Public Member Functions

 ExodusII_IO_Helper (const ParallelObject &parent, bool v=false, bool run_only_on_proc0=true, bool single_precision=false)
 
virtual ~ExodusII_IO_Helper ()
 
const char * get_elem_type () const
 
void open (const char *filename, bool read_only)
 
void read_header ()
 
void read_qa_records ()
 
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, std::map< dof_id_type, Real > &elem_var_value_map)
 
virtual void create (std::string filename)
 
virtual void initialize (std::string title, const MeshBase &mesh, bool use_discontinuous=false)
 
virtual void write_nodal_coordinates (const MeshBase &mesh, bool use_discontinuous=false)
 
virtual void write_elements (const MeshBase &mesh, bool use_discontinuous=false)
 
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< Real > &values, int timestep)
 
void write_nodal_values (int var_id, const std::vector< Real > &values, int timestep)
 
void write_information_records (const std::vector< std::string > &records)
 
void write_global_values (const std::vector< Real > &values, int timestep)
 
void use_mesh_dimension_instead_of_spatial_dimension (bool val)
 
void write_as_dimension (unsigned dim)
 
void set_coordinate_offset (Point p)
 
std::vector< std::string > get_complex_names (const std::vector< std::string > &names) const
 
void message (const std::string &msg)
 
void message (const std::string &msg, int i)
 
void read_var_names (ExodusVarType type)
 
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
 
unsigned _write_as_dimension
 
Point _coordinate_offset
 
bool _single_precision
 
const Parallel::Communicator_communicator
 

Private Member Functions

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
Date
2002

Definition at line 83 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 570 of file exodusII_io_helper.h.

Constructor & Destructor Documentation

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

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 265 of file exodusII_io_helper.C.

References elem_type, and title.

268  :
269  ParallelObject(parent),
270  ex_id(0),
271  ex_err(0),
272  num_dim(0),
273  num_global_vars(0),
274  num_nodes(0),
275  num_elem(0),
276  num_elem_blk(0),
277  num_node_sets(0),
278  num_side_sets(0),
281  num_attr(0),
283  num_time_steps(0),
284  num_nodal_vars(0),
285  num_elem_vars(0),
286  verbose(v),
287  opened_for_writing(false),
288  opened_for_reading(false),
289  _run_only_on_proc0(run_only_on_proc0),
290  _elem_vars_initialized(false),
295  _single_precision(single_precision)
296 {
297  title.resize(MAX_LINE_LENGTH+1);
298  elem_type.resize(MAX_STR_LENGTH);
299 }
ParallelObject(const Parallel::Communicator &comm_in)
libMesh::ExodusII_IO_Helper::~ExodusII_IO_Helper ( )
virtual

Destructor.

Definition at line 303 of file exodusII_io_helper.C.

304 {
305 }

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 1761 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().

1764 {
1765  // There may already be global variables in the file (for example,
1766  // if we're appending) and in that case, we
1767  // 1.) Cannot initialize them again.
1768  // 2.) Should check to be sure that the global variable names are the same.
1769 
1770  // Fills up names_from_file for us
1771  this->read_var_names(type);
1772 
1773  // Both the names of the global variables and their order must match
1774  if (names_from_file != names)
1775  {
1776  libMesh::err << "Error! The Exodus file already contains the variables:" << std::endl;
1777  for (std::size_t i=0; i<names_from_file.size(); ++i)
1778  libMesh::out << names_from_file[i] << std::endl;
1779 
1780  libMesh::err << "And you asked to write:" << std::endl;
1781  for (std::size_t i=0; i<names.size(); ++i)
1782  libMesh::out << names[i] << std::endl;
1783 
1784  libmesh_error_msg("Cannot overwrite existing variables in Exodus II file.");
1785  }
1786 }
void read_var_names(ExodusVarType type)
OStreamProxy err(std::cerr)
OStreamProxy out(std::cout)
void libMesh::ExodusII_IO_Helper::close ( )

Closes the ExodusII mesh file.

Definition at line 782 of file exodusII_io_helper.C.

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

Referenced by libMesh::Nemesis_IO_Helper::~Nemesis_IO_Helper().

783 {
784  // Always call close on processor 0.
785  // If we're running on multiple processors, i.e. as one of several Nemesis files,
786  // we call close on all processors...
787  if ((this->processor_id() == 0) || (!_run_only_on_proc0))
788  {
789  // Don't close the file if it was never opened, this raises an Exodus error
791  {
793  EX_CHECK_ERR(ex_err, "Error closing Exodus file.");
794  message("Exodus file closed successfully.");
795  }
796  }
797 }
void message(const std::string &msg)
EXODUS_EXPORT int ex_close(int exoid)
processor_id_type processor_id() const
const Parallel::Communicator& libMesh::ParallelObject::comm ( ) const
inlineinherited
Returns
a reference to the Parallel::Communicator object used by this mesh.

Definition at line 87 of file parallel_object.h.

References libMesh::ParallelObject::_communicator.

Referenced by libMesh::__libmesh_petsc_diff_solver_monitor(), libMesh::__libmesh_petsc_diff_solver_residual(), libMesh::__libmesh_petsc_snes_jacobian(), libMesh::__libmesh_petsc_snes_postcheck(), libMesh::__libmesh_petsc_snes_residual(), libMesh::__libmesh_tao_equality_constraints(), libMesh::__libmesh_tao_equality_constraints_jacobian(), libMesh::__libmesh_tao_gradient(), libMesh::__libmesh_tao_hessian(), libMesh::__libmesh_tao_inequality_constraints(), libMesh::__libmesh_tao_inequality_constraints_jacobian(), libMesh::__libmesh_tao_objective(), libMesh::MeshRefinement::_coarsen_elements(), libMesh::ExactSolution::_compute_error(), libMesh::ParmetisPartitioner::_do_repartition(), libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::BoundaryInfo::_find_id_maps(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_get_diagonal(), libMesh::SlepcEigenSolver< T >::_petsc_shell_matrix_get_diagonal(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_mult(), libMesh::SlepcEigenSolver< T >::_petsc_shell_matrix_mult(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_mult_add(), libMesh::EquationSystems::_read_impl(), libMesh::MeshRefinement::_refine_elements(), libMesh::MeshRefinement::_smooth_flags(), libMesh::ImplicitSystem::add_matrix(), libMesh::System::add_vector(), libMesh::EigenSparseLinearSolver< T >::adjoint_solve(), libMesh::UnstructuredMesh::all_second_order(), libMesh::MeshTools::Modification::all_tri(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::FEMSystem::assemble_qoi(), libMesh::MeshCommunication::assign_global_indices(), libMesh::ParmetisPartitioner::assign_partitioning(), libMesh::DofMap::attach_matrix(), libMesh::Parallel::BinSorter< KeyType, IdxType >::binsort(), libMesh::Parallel::Sort< KeyType, IdxType >::binsort(), libMesh::MeshTools::bounding_box(), libMesh::MeshCommunication::broadcast(), libMesh::SparseMatrix< T >::build(), libMesh::MeshTools::Generation::build_extrusion(), libMesh::Parallel::Histogram< KeyType, IdxType >::build_histogram(), libMesh::PetscNonlinearSolver< T >::build_mat_null_space(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::EquationSystems::build_parallel_solution_vector(), libMesh::MeshBase::cache_elem_dims(), libMesh::System::calculate_norm(), libMesh::DofMap::check_dirichlet_bcid_consistency(), libMesh::DistributedVector< T >::clone(), libMesh::EigenSparseVector< T >::clone(), libMesh::LaspackVector< T >::clone(), libMesh::EpetraVector< T >::clone(), libMesh::PetscVector< T >::clone(), libMesh::EpetraVector< T >::close(), libMesh::Parallel::Sort< KeyType, IdxType >::communicate_bins(), 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::MeshRefinement::create_parent_error_vector(), libMesh::MeshCommunication::delete_remote_elements(), libMesh::DofMap::distribute_dofs(), DMlibMeshFunction(), DMlibMeshJacobian(), DMlibMeshSetSystem_libMesh(), DMVariableBounds_libMesh(), libMesh::MeshRefinement::eliminate_unrefined_patches(), libMesh::EpetraVector< T >::EpetraVector(), libMesh::WeightedPatchRecoveryErrorEstimator::estimate_error(), libMesh::PatchRecoveryErrorEstimator::estimate_error(), libMesh::JumpErrorEstimator::estimate_error(), libMesh::AdjointRefinementEstimator::estimate_error(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::MeshRefinement::flag_elements_by_elem_fraction(), libMesh::MeshRefinement::flag_elements_by_error_fraction(), libMesh::MeshRefinement::flag_elements_by_nelem_target(), libMesh::MeshCommunication::gather(), libMesh::MeshCommunication::gather_neighboring_elements(), libMesh::CondensedEigenSystem::get_eigenpair(), libMesh::DofMap::get_info(), libMesh::ImplicitSystem::get_linear_solver(), libMesh::EquationSystems::get_solution(), libMesh::LocationMap< T >::init(), libMesh::PetscDiffSolver::init(), libMesh::TimeSolver::init(), libMesh::TopologyMap::init(), libMesh::TaoOptimizationSolver< T >::init(), libMesh::PetscNonlinearSolver< T >::init(), libMesh::DistributedVector< T >::init(), libMesh::EpetraVector< T >::init(), libMesh::PetscVector< T >::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::EigenSystem::init_data(), libMesh::EigenSystem::init_matrices(), libMesh::ParmetisPartitioner::initialize(), libMesh::OptimizationSystem::initialize_equality_constraints_storage(), libMesh::OptimizationSystem::initialize_inequality_constraints_storage(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Elem >(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_topology_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_valid_boundary_ids(), libMesh::MeshTools::libmesh_assert_valid_dof_ids(), libMesh::MeshTools::libmesh_assert_valid_neighbors(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_flags(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_object_ids(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_p_levels(), libMesh::MeshTools::libmesh_assert_valid_refinement_flags(), libMesh::MeshTools::libmesh_assert_valid_unique_ids(), libMesh::MeshRefinement::limit_level_mismatch_at_edge(), libMesh::MeshRefinement::limit_level_mismatch_at_node(), libMesh::MeshRefinement::limit_overrefined_boundary(), libMesh::MeshRefinement::limit_underrefined_boundary(), libMesh::MeshRefinement::make_coarsening_compatible(), libMesh::MeshCommunication::make_elems_parallel_consistent(), libMesh::MeshRefinement::make_flags_parallel_consistent(), libMesh::MeshCommunication::make_new_node_proc_ids_parallel_consistent(), libMesh::MeshCommunication::make_new_nodes_parallel_consistent(), libMesh::MeshCommunication::make_node_ids_parallel_consistent(), libMesh::MeshCommunication::make_node_proc_ids_parallel_consistent(), libMesh::MeshCommunication::make_node_unique_ids_parallel_consistent(), libMesh::MeshCommunication::make_nodes_parallel_consistent(), libMesh::MeshCommunication::make_p_levels_parallel_consistent(), libMesh::MeshRefinement::make_refinement_compatible(), libMesh::DistributedVector< T >::max(), libMesh::FEMSystem::mesh_position_set(), libMesh::MeshSerializer::MeshSerializer(), libMesh::DistributedVector< T >::min(), libMesh::DistributedMesh::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::BoundaryInfo::n_shellface_conds(), libMesh::DistributedMesh::parallel_max_elem_id(), libMesh::DistributedMesh::parallel_max_node_id(), libMesh::ReplicatedMesh::parallel_max_unique_id(), libMesh::DistributedMesh::parallel_max_unique_id(), libMesh::DistributedMesh::parallel_n_elem(), libMesh::DistributedMesh::parallel_n_nodes(), libMesh::SparsityPattern::Build::parallel_sync(), libMesh::MeshTools::paranoid_n_levels(), libMesh::Partitioner::partition(), libMesh::MetisPartitioner::partition_range(), libMesh::Partitioner::partition_unpartitioned_elements(), libMesh::petsc_auto_fieldsplit(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::MeshBase::prepare_for_use(), libMesh::SparseMatrix< T >::print(), libMesh::MeshTools::processor_bounding_box(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), libMesh::XdrIO::read_serialized_bc_names(), libMesh::XdrIO::read_serialized_bcs_helper(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::MeshBase::recalculate_n_partitions(), libMesh::MeshCommunication::redistribute(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::MeshCommunication::send_coarse_ghosts(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::Partitioner::set_parent_processor_ids(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::Parallel::Sort< KeyType, IdxType >::sort(), libMesh::MeshTools::subdomain_bounding_box(), libMesh::MeshBase::subdomain_ids(), libMesh::BoundaryInfo::sync(), libMesh::Parallel::sync_element_data_by_parent_id(), libMesh::Parallel::sync_node_data_by_element_id(), libMesh::MeshRefinement::test_level_one(), libMesh::MeshRefinement::test_unflagged(), libMesh::MeshTools::total_weight(), libMesh::NameBasedIO::write(), libMesh::XdrIO::write(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), libMesh::XdrIO::write_serialized_nodesets(), libMesh::DistributedVector< T >::zero_clone(), libMesh::LaspackVector< T >::zero_clone(), libMesh::EigenSparseVector< T >::zero_clone(), libMesh::EpetraVector< T >::zero_clone(), and libMesh::PetscVector< T >::zero_clone().

88  { return _communicator; }
const Parallel::Communicator & _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 1066 of file exodusII_io_helper.C.

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

1067 {
1068  // If we're processor 0, always create the file.
1069  // If we running on all procs, e.g. as one of several Nemesis files, also
1070  // call create there.
1071  if ((this->processor_id() == 0) || (!_run_only_on_proc0))
1072  {
1073  int
1074  comp_ws = 0,
1075  io_ws = 0;
1076 
1077  if(_single_precision)
1078  {
1079  comp_ws = cast_int<int>(sizeof(float));
1080  io_ws = cast_int<int>(sizeof(float));
1081  }
1082  // Fall back on double precision when necessary since ExodusII
1083  // doesn't seem to support long double
1084  else
1085  {
1086  comp_ws = cast_int<int>
1087  (std::min(sizeof(Real), sizeof(double)));
1088  io_ws = cast_int<int>
1089  (std::min(sizeof(Real), sizeof(double)));
1090  }
1091 
1092  ex_id = exII::ex_create(filename.c_str(), EX_CLOBBER, &comp_ws, &io_ws);
1093 
1094  EX_CHECK_ERR(ex_id, "Error creating ExodusII mesh file.");
1095 
1096  if (verbose)
1097  libMesh::out << "File created successfully." << std::endl;
1098  }
1099 
1100  opened_for_writing = true;
1101  current_filename = filename;
1102 }
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
OStreamProxy out(std::cout)
long double min(long double a, double b)
processor_id_type processor_id() const
int libMesh::ExodusII_IO_Helper::get_block_id ( int  index)

Get the block number for the given block index.

Definition at line 531 of file exodusII_io_helper.C.

References block_ids.

Referenced by write_element_values().

532 {
533  libmesh_assert_less (static_cast<unsigned int>(index), block_ids.size());
534 
535  return block_ids[index];
536 }
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 540 of file exodusII_io_helper.C.

References block_ids, and id_to_block_names.

541 {
542  libmesh_assert_less (static_cast<unsigned int>(index), block_ids.size());
543 
544  return id_to_block_names[block_ids[index]];
545 }
std::map< int, std::string > id_to_block_names
std::vector< std::string > libMesh::ExodusII_IO_Helper::get_complex_names ( const std::vector< std::string > &  names) const

Returns a vector with three copies of each element in the provided name vector, starting with r_, i_ and a_ respectively.

Definition at line 1982 of file exodusII_io_helper.C.

1983 {
1984  std::vector<std::string>::const_iterator names_it = names.begin();
1985  std::vector<std::string>::const_iterator names_end = names.end();
1986 
1987  std::vector<std::string> complex_names;
1988 
1989  // This will loop over all names and create new "complex" names
1990  // (i.e. names that start with r_, i_ or a_
1991  for (; names_it != names_end; ++names_it)
1992  {
1993  std::stringstream name_real, name_imag, name_abs;
1994  name_real << "r_" << *names_it;
1995  name_imag << "i_" << *names_it;
1996  name_abs << "a_" << *names_it;
1997 
1998  complex_names.push_back(name_real.str());
1999  complex_names.push_back(name_imag.str());
2000  complex_names.push_back(name_abs.str());
2001  }
2002 
2003  return complex_names;
2004 }
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 309 of file exodusII_io_helper.C.

References elem_type.

310 {
311  return &elem_type[0];
312 }
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 567 of file exodusII_io_helper.C.

References nodeset_ids.

568 {
569  libmesh_assert_less (static_cast<unsigned int>(index), nodeset_ids.size());
570 
571  return nodeset_ids[index];
572 }
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 576 of file exodusII_io_helper.C.

References id_to_ns_names, and nodeset_ids.

577 {
578  libmesh_assert_less (static_cast<unsigned int>(index), nodeset_ids.size());
579 
580  return id_to_ns_names[nodeset_ids[index]];
581 }
std::map< int, std::string > id_to_ns_names
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 549 of file exodusII_io_helper.C.

References ss_ids.

550 {
551  libmesh_assert_less (static_cast<unsigned int>(index), ss_ids.size());
552 
553  return ss_ids[index];
554 }
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 558 of file exodusII_io_helper.C.

References id_to_ss_names, and ss_ids.

559 {
560  libmesh_assert_less (static_cast<unsigned int>(index), ss_ids.size());
561 
562  return id_to_ss_names[ss_ids[index]];
563 }
std::map< int, std::string > id_to_ss_names
void libMesh::ExodusII_IO_Helper::initialize ( std::string  title,
const MeshBase mesh,
bool  use_discontinuous = false 
)
virtual

Initializes the Exodus file.

Reimplemented in libMesh::Nemesis_IO_Helper.

Definition at line 1106 of file exodusII_io_helper.C.

References _run_only_on_proc0, _use_mesh_dimension_instead_of_spatial_dimension, _write_as_dimension, libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::BoundaryInfo::build_node_boundary_ids(), libMesh::BoundaryInfo::build_shellface_boundary_ids(), libMesh::BoundaryInfo::build_side_boundary_ids(), end, libMesh::err, ex_err, ex_id, ex_put_init(), libMesh::MeshBase::get_boundary_info(), libMesh::DofObject::id(), libMesh::MeshBase::mesh_dimension(), libMesh::MeshBase::n_active_elem(), libMesh::MeshBase::n_elem(), libMesh::MeshBase::nodes_begin(), libMesh::MeshBase::nodes_end(), 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().

1107 {
1108  // n_active_elem() is a parallel_only function
1109  unsigned int n_active_elem = mesh.n_active_elem();
1110 
1111  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1112  return;
1113 
1114  // If _write_as_dimension is nonzero, use it to set num_dim in the Exodus file.
1115  if (_write_as_dimension)
1118  num_dim = mesh.mesh_dimension();
1119  else
1120  num_dim = mesh.spatial_dimension();
1121 
1122  num_elem = mesh.n_elem();
1123 
1124  if (!use_discontinuous)
1125  {
1126  // Don't rely on mesh.n_nodes() here. If nodes have been
1127  // deleted recently, it will be incorrect.
1128  num_nodes = cast_int<int>(std::distance(mesh.nodes_begin(),
1129  mesh.nodes_end()));
1130  }
1131  else
1132  {
1135  for (; it!=end; ++it)
1136  num_nodes += (*it)->n_nodes();
1137  }
1138 
1139  std::vector<boundary_id_type> unique_side_boundaries;
1140  std::vector<boundary_id_type> unique_node_boundaries;
1141 
1142  mesh.get_boundary_info().build_side_boundary_ids(unique_side_boundaries);
1143  {
1144  // Add shell face boundaries to the list of side boundaries, since ExodusII
1145  // treats these the same way.
1146  std::vector<boundary_id_type> shellface_boundaries;
1147  mesh.get_boundary_info().build_shellface_boundary_ids(shellface_boundaries);
1148  for (std::size_t i=0; i<shellface_boundaries.size(); i++)
1149  unique_side_boundaries.push_back(shellface_boundaries[i]);
1150  }
1151  mesh.get_boundary_info().build_node_boundary_ids(unique_node_boundaries);
1152 
1153  num_side_sets = cast_int<int>(unique_side_boundaries.size());
1154  num_node_sets = cast_int<int>(unique_node_boundaries.size());
1155 
1156  //loop through element and map between block and element vector
1157  std::map<subdomain_id_type, std::vector<unsigned int> > subdomain_map;
1158 
1161  for (; it!=end; ++it)
1162  {
1163  const Elem * elem = *it;
1164  subdomain_id_type cur_subdomain = elem->subdomain_id();
1165 
1166  subdomain_map[cur_subdomain].push_back(elem->id());
1167  }
1168  num_elem_blk = cast_int<int>(subdomain_map.size());
1169 
1170  if (str_title.size() > MAX_LINE_LENGTH)
1171  {
1172  libMesh::err << "Warning, Exodus files cannot have titles longer than "
1173  << MAX_LINE_LENGTH
1174  << " characters. Your title will be truncated."
1175  << std::endl;
1176  str_title.resize(MAX_LINE_LENGTH);
1177  }
1178 
1180  str_title.c_str(),
1181  num_dim,
1182  num_nodes,
1183  n_active_elem,
1184  num_elem_blk,
1185  num_node_sets,
1186  num_side_sets);
1187 
1188  EX_CHECK_ERR(ex_err, "Error initializing new Exodus file.");
1189 }
const BoundaryInfo & get_boundary_info() const
Definition: mesh_base.h:111
virtual dof_id_type n_active_elem() const =0
The base class for all geometric element types.
Definition: elem.h:86
IterBase * end
EXODUS_EXPORT int ex_put_init(int exoid, const char *title, int64_t num_dim, int64_t num_nodes, int64_t num_elem, int64_t num_elem_blk, int64_t num_node_sets, int64_t num_side_sets)
virtual node_iterator nodes_begin()=0
void build_side_boundary_ids(std::vector< boundary_id_type > &b_ids) const
unsigned int spatial_dimension() const
Definition: mesh_base.C:157
OStreamProxy err(std::cerr)
subdomain_id_type subdomain_id() const
Definition: elem.h:1733
virtual element_iterator active_elements_begin()=0
virtual element_iterator active_elements_end()=0
virtual node_iterator nodes_end()=0
unsigned int mesh_dimension() const
Definition: mesh_base.C:148
void build_shellface_boundary_ids(std::vector< boundary_id_type > &b_ids) const
dof_id_type id() const
Definition: dof_object.h:624
virtual dof_id_type n_elem() const =0
void build_node_boundary_ids(std::vector< boundary_id_type > &b_ids) const
processor_id_type processor_id() const
void libMesh::ExodusII_IO_Helper::initialize_element_variables ( std::vector< std::string >  names)

Sets up the nodal variables

Definition at line 1660 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().

1661 {
1662  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1663  return;
1664 
1665  // Quick return if there are no element variables to write
1666  if (names.size() == 0)
1667  return;
1668 
1669  // Quick return if we have already called this function
1671  return;
1672 
1673  // Be sure that variables in the file match what we are asking for
1674  if (num_elem_vars > 0)
1675  {
1676  this->check_existing_vars(ELEMENTAL, names, this->elem_var_names);
1677  return;
1678  }
1679 
1680  // Set the flag so we can skip this stuff on subsequent calls to
1681  // initialize_element_variables()
1682  _elem_vars_initialized = true;
1683 
1684  this->write_var_names(ELEMENTAL, names);
1685 
1686  // Form the element variable truth table and send to Exodus.
1687  // This tells which variables are written to which blocks,
1688  // and can dramatically speed up writing element variables
1689  //
1690  // We really should initialize all entries in the truth table to 0
1691  // and then loop over all subdomains, setting their entries to 1
1692  // if a given variable exists on that subdomain. However,
1693  // we don't have that information, and the element variables
1694  // passed to us are padded with zeroes for the blocks where
1695  // they aren't defined. To be consistent with that, fill
1696  // the truth table with ones.
1697  std::vector<int> truth_tab(num_elem_blk*num_elem_vars, 1);
1699  num_elem_blk,
1700  num_elem_vars,
1701  &truth_tab[0]);
1702  EX_CHECK_ERR(ex_err, "Error writing element truth table.");
1703 }
std::vector< std::string > elem_var_names
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)
EXODUS_EXPORT int ex_put_elem_var_tab(int exoid, int num_elem_blk, int num_elem_var, int *elem_var_tab)
processor_id_type processor_id() const
void libMesh::ExodusII_IO_Helper::initialize_global_variables ( std::vector< std::string >  names)

Sets up the global variables

Definition at line 1735 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().

1736 {
1737  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1738  return;
1739 
1740  // Quick return if there are no global variables to write
1741  if (names.size() == 0)
1742  return;
1743 
1745  return;
1746 
1747  // Be sure that variables in the file match what we are asking for
1748  if (num_global_vars > 0)
1749  {
1750  this->check_existing_vars(GLOBAL, names, this->global_var_names);
1751  return;
1752  }
1753 
1754  _global_vars_initialized = true;
1755 
1756  this->write_var_names(GLOBAL, names);
1757 }
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)
std::vector< std::string > global_var_names
processor_id_type processor_id() const
void libMesh::ExodusII_IO_Helper::initialize_nodal_variables ( std::vector< std::string >  names)

Sets up the nodal variables

Definition at line 1707 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().

1708 {
1709  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1710  return;
1711 
1712  // Quick return if there are no nodal variables to write
1713  if (names.size() == 0)
1714  return;
1715 
1716  // Quick return if we have already called this function
1718  return;
1719 
1720  // Be sure that variables in the file match what we are asking for
1721  if (num_nodal_vars > 0)
1722  {
1723  this->check_existing_vars(NODAL, names, this->nodal_var_names);
1724  return;
1725  }
1726 
1727  // Set the flag so we can skip the rest of this function on subsequent calls.
1728  _nodal_vars_initialized = true;
1729 
1730  this->write_var_names(NODAL, names);
1731 }
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)
std::vector< std::string > nodal_var_names
processor_id_type processor_id() const
int libMesh::ExodusII_IO_Helper::inquire ( int  req_info,
std::string  error_msg = "" 
)

Generic inquiry, returns the value

Definition at line 801 of file exodusII_io_helper.C.

References ex_err, ex_id, and ex_inquire().

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

802 {
803  int ret_int = 0;
804  char ret_char = 0;
805  float ret_float = 0.;
806 
808  req_info_in,
809  &ret_int,
810  &ret_float,
811  &ret_char);
812 
813  EX_CHECK_ERR(ex_err, error_msg);
814 
815  return ret_int;
816 }
EXODUS_EXPORT int ex_inquire(int exoid, int inquiry, void_int *, float *, char *)
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 316 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().

317 {
318  if (verbose) libMesh::out << msg << std::endl;
319 }
OStreamProxy out(std::cout)
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 323 of file exodusII_io_helper.C.

References libMesh::out, and verbose.

324 {
325  if (verbose) libMesh::out << msg << i << "." << std::endl;
326 }
OStreamProxy out(std::cout)
processor_id_type libMesh::ParallelObject::n_processors ( ) const
inlineinherited
Returns
the number of processors in the group.

Definition at line 93 of file parallel_object.h.

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

Referenced by libMesh::ParmetisPartitioner::_do_repartition(), libMesh::BoundaryInfo::_find_id_maps(), libMesh::DistributedMesh::add_elem(), libMesh::DistributedMesh::add_node(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::FEMSystem::assembly(), libMesh::ParmetisPartitioner::assign_partitioning(), libMesh::AztecLinearSolver< T >::AztecLinearSolver(), libMesh::MeshCommunication::broadcast(), libMesh::BoundaryInfo::build_node_list_from_side_list(), libMesh::DistributedMesh::clear(), libMesh::Nemesis_IO_Helper::compute_border_node_ids(), libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::MeshTools::correct_node_proc_ids(), libMesh::UnstructuredMesh::create_pid_mesh(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::DistributedMesh::DistributedMesh(), libMesh::EnsightIO::EnsightIO(), libMesh::MeshCommunication::gather(), libMesh::MeshCommunication::gather_neighboring_elements(), libMesh::MeshBase::get_info(), libMesh::EquationSystems::get_solution(), libMesh::DistributedVector< T >::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::ParmetisPartitioner::initialize(), libMesh::Nemesis_IO_Helper::initialize(), libMesh::DistributedMesh::insert_elem(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Elem >(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_topology_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_valid_boundary_ids(), libMesh::MeshTools::libmesh_assert_valid_dof_ids(), libMesh::MeshTools::libmesh_assert_valid_neighbors(), libMesh::MeshTools::libmesh_assert_valid_refinement_flags(), libMesh::DofMap::local_variable_indices(), libMesh::MeshBase::n_active_elem_on_proc(), libMesh::MeshBase::n_elem_on_proc(), libMesh::MeshBase::n_nodes_on_proc(), libMesh::SparsityPattern::Build::parallel_sync(), 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::SparseMatrix< T >::print(), libMesh::MeshTools::processor_bounding_box(), libMesh::NameBasedIO::read(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), libMesh::MeshCommunication::redistribute(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::Partitioner::repartition(), libMesh::MeshCommunication::send_coarse_ghosts(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::Parallel::Sort< KeyType, IdxType >::sort(), libMesh::DistributedMesh::update_parallel_id_counts(), libMesh::GMVIO::write_binary(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), and libMesh::XdrIO::write_serialized_nodesets().

94  { return cast_int<processor_id_type>(_communicator.size()); }
unsigned int size() const
Definition: parallel.h:679
const Parallel::Communicator & _communicator
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 330 of file exodusII_io_helper.C.

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

331 {
332  // Version of Exodus you are using
333  float ex_version = 0.;
334 
335  // Word size in bytes of the floating point variables used in the
336  // application program. Exodus only supports 4-byte and 8-byte
337  // floats.
338  int comp_ws = std::min(sizeof(Real), std::size_t(8));
339 
340  // Word size in bytes of the floating point data as they are stored
341  // in the ExodusII file. "If this argument is 0, the word size of the
342  // floating point data already stored in the file is returned"
343  int io_ws = 0;
344 
345  ex_id = exII::ex_open(filename,
346  read_only ? EX_READ : EX_WRITE,
347  &comp_ws,
348  &io_ws,
349  &ex_version);
350 
351  std::string err_msg = std::string("Error opening ExodusII mesh file: ") + std::string(filename);
352  EX_CHECK_ERR(ex_id, err_msg);
353  if (verbose) libMesh::out << "File opened successfully." << std::endl;
354 
355  if (read_only)
356  opened_for_reading = true;
357  else
358  opened_for_writing = true;
359 
360  current_filename = std::string(filename);
361 }
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
OStreamProxy out(std::cout)
long double min(long double a, double b)
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 445 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.

446 {
447  if (verbose)
448  libMesh::out << "Title: \t" << &title[0] << std::endl
449  << "Mesh Dimension: \t" << num_dim << std::endl
450  << "Number of Nodes: \t" << num_nodes << std::endl
451  << "Number of elements: \t" << num_elem << std::endl
452  << "Number of elt blocks: \t" << num_elem_blk << std::endl
453  << "Number of node sets: \t" << num_node_sets << std::endl
454  << "Number of side sets: \t" << num_side_sets << std::endl;
455 }
OStreamProxy out(std::cout)
void libMesh::ExodusII_IO_Helper::print_nodes ( std::ostream &  out = libMesh::out)

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

Definition at line 496 of file exodusII_io_helper.C.

References num_nodes, x, y, and z.

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

Definition at line 99 of file parallel_object.h.

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

Referenced by libMesh::BoundaryInfo::_find_id_maps(), libMesh::EquationSystems::_read_impl(), libMesh::DistributedMesh::add_elem(), libMesh::BoundaryInfo::add_elements(), libMesh::DofMap::add_neighbors_to_send_list(), libMesh::DistributedMesh::add_node(), libMesh::MeshRefinement::add_node(), libMesh::MeshTools::Modification::all_tri(), libMesh::FEMSystem::assembly(), libMesh::ParmetisPartitioner::assign_partitioning(), libMesh::MeshCommunication::broadcast(), 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::BoundaryInfo::build_node_list_from_side_list(), libMesh::DofMap::build_sparsity(), libMesh::DistributedMesh::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(), libMesh::MeshTools::correct_node_proc_ids(), create(), libMesh::DistributedMesh::delete_elem(), libMesh::DistributedMesh::delete_node(), libMesh::MeshCommunication::delete_remote_elements(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::DistributedMesh::DistributedMesh(), libMesh::EnsightIO::EnsightIO(), libMesh::MeshFunction::find_element(), libMesh::MeshFunction::find_elements(), libMesh::UnstructuredMesh::find_neighbors(), libMesh::MeshCommunication::gather(), libMesh::MeshCommunication::gather_neighboring_elements(), 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::DofMap::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::EquationSystems::get_solution(), libMesh::Nemesis_IO_Helper::get_ss_param_global(), libMesh::SparsityPattern::Build::handle_vi_vj(), libMesh::DistributedVector< T >::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::ParmetisPartitioner::initialize(), initialize(), initialize_element_variables(), initialize_global_variables(), initialize_nodal_variables(), libMesh::DistributedMesh::insert_elem(), libMesh::DofMap::is_evaluable(), libMesh::SparsityPattern::Build::join(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Elem >(), libMesh::MeshTools::libmesh_assert_parallel_consistent_procids< Node >(), libMesh::MeshTools::libmesh_assert_valid_neighbors(), libMesh::DistributedMesh::libmesh_assert_valid_parallel_object_ids(), libMesh::DofMap::local_variable_indices(), libMesh::MeshRefinement::make_coarsening_compatible(), libMesh::MeshBase::n_active_local_elem(), libMesh::BoundaryInfo::n_boundary_conds(), libMesh::BoundaryInfo::n_edge_conds(), libMesh::System::n_local_dofs(), libMesh::MeshBase::n_local_elem(), libMesh::MeshBase::n_local_nodes(), libMesh::BoundaryInfo::n_nodeset_conds(), libMesh::BoundaryInfo::n_shellface_conds(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::SparsityPattern::Build::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::SparsityPattern::Build::parallel_sync(), libMesh::MetisPartitioner::partition_range(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::SparseMatrix< T >::print(), libMesh::NumericVector< T >::print_global(), 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::NameBasedIO::read(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), read_elem_num_map(), read_node_num_map(), libMesh::XdrIO::read_serialized_bc_names(), libMesh::XdrIO::read_serialized_bcs_helper(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::MeshCommunication::redistribute(), libMesh::DistributedMesh::renumber_dof_objects(), libMesh::MeshCommunication::send_coarse_ghosts(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::MeshTools::total_weight(), libMesh::Parallel::Packing< Node * >::unpack(), libMesh::Parallel::Packing< Elem * >::unpack(), libMesh::DistributedMesh::update_parallel_id_counts(), libMesh::MeshTools::weight(), libMesh::NameBasedIO::write(), libMesh::XdrIO::write(), libMesh::EquationSystems::write(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::ExodusII_IO::write_element_data(), write_element_values(), write_elements(), libMesh::ExodusII_IO::write_global_data(), write_global_values(), libMesh::ExodusII_IO::write_information_records(), write_information_records(), write_nodal_coordinates(), libMesh::UCDIO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data_discontinuous(), write_nodal_values(), libMesh::Nemesis_IO_Helper::write_nodesets(), write_nodesets(), libMesh::XdrIO::write_serialized_bc_names(), libMesh::XdrIO::write_serialized_bcs_helper(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), libMesh::XdrIO::write_serialized_nodesets(), libMesh::XdrIO::write_serialized_subdomain_names(), libMesh::Nemesis_IO_Helper::write_sidesets(), write_sidesets(), libMesh::ExodusII_IO::write_timestep(), and write_timestep().

100  { return cast_int<processor_id_type>(_communicator.rank()); }
const Parallel::Communicator & _communicator
unsigned int rank() const
Definition: parallel.h:677
void libMesh::ExodusII_IO_Helper::read_block_info ( )

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

Definition at line 504 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, libmesh_nullptr, message(), and num_elem_blk.

505 {
506  block_ids.resize(num_elem_blk);
507  // Get all element block IDs.
509  block_ids.empty() ? libmesh_nullptr : &block_ids[0]);
510  // Usually, there is only one
511  // block since there is only
512  // one type of element.
513  // However, there could be more.
514 
515  EX_CHECK_ERR(ex_err, "Error getting block IDs.");
516  message("All block IDs retrieved successfully.");
517 
518  char name_buffer[MAX_STR_LENGTH+1];
519  for (int i=0; i<num_elem_blk; ++i)
520  {
522  block_ids[i], name_buffer);
523  EX_CHECK_ERR(ex_err, "Error getting block name.");
524  id_to_block_names[block_ids[i]] = name_buffer;
525  }
526  message("All block names retrieved successfully.");
527 }
EXODUS_EXPORT int ex_get_name(int exoid, ex_entity_type obj_type, ex_entity_id entity_id, char *name)
const class libmesh_nullptr_t libmesh_nullptr
EXODUS_EXPORT int ex_get_elem_blk_ids(int exoid, void_int *ids)
void message(const std::string &msg)
std::map< int, std::string > id_to_block_names
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 586 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.

587 {
588  libmesh_assert_less (static_cast<unsigned int>(block), block_ids.size());
589 
591  block_ids[block],
592  &elem_type[0],
595  &num_attr);
596  if (verbose)
597  libMesh::out << "Reading a block of " << num_elem_this_blk
598  << " " << &elem_type[0] << "(s)"
599  << " having " << num_nodes_per_elem
600  << " nodes per element." << std::endl;
601 
602  EX_CHECK_ERR(ex_err, "Error getting block info.");
603  message("Info retrieved successfully for block: ", block);
604 
605 
606 
607  // Read in the connectivity of the elements of this block,
608  // watching out for the case where we actually have no
609  // elements in this block (possible with parallel files)
611 
612  if (!connect.empty())
613  {
614  ex_err = exII::ex_get_elem_conn(ex_id,
615  block_ids[block],
616  &connect[0]);
617 
618  EX_CHECK_ERR(ex_err, "Error reading block connectivity.");
619  message("Connectivity retrieved successfully for block: ", block);
620  }
621 }
EXODUS_EXPORT int ex_get_elem_block(int exoid, ex_entity_id elem_blk_id, char *elem_type, void_int *num_elem_this_blk, void_int *num_nodes_per_elem, void_int *num_attr)
void message(const std::string &msg)
EXODUS_EXPORT int ex_get_elem_conn(int exoid, ex_entity_id elem_blk_id, void_int *connect)
OStreamProxy out(std::cout)
void libMesh::ExodusII_IO_Helper::read_elem_num_map ( )

Reads the optional node_num_map from the ExodusII mesh file.

Definition at line 626 of file exodusII_io_helper.C.

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

627 {
628  elem_num_map.resize(num_elem);
629 
631  elem_num_map.empty() ? libmesh_nullptr : &elem_num_map[0]);
632 
633  EX_CHECK_ERR(ex_err, "Error retrieving element number map.");
634  message("Element numbering map retrieved successfully.");
635 
636 
637  if (verbose)
638  {
639  libMesh::out << "[" << this->processor_id() << "] elem_num_map[i] = ";
640  for (unsigned int i=0; i<static_cast<unsigned int>(std::min(10, num_elem-1)); ++i)
641  libMesh::out << elem_num_map[i] << ", ";
642  libMesh::out << "... " << elem_num_map.back() << std::endl;
643  }
644 }
const class libmesh_nullptr_t libmesh_nullptr
EXODUS_EXPORT int ex_get_elem_num_map(int exoid, void_int *elem_map)
void message(const std::string &msg)
OStreamProxy out(std::cout)
long double min(long double a, double b)
processor_id_type processor_id() const
void libMesh::ExodusII_IO_Helper::read_elemental_var_values ( std::string  elemental_var_name,
int  time_step,
std::map< dof_id_type, Real > &  elem_var_value_map 
)

Reads elemental values for the variable 'elemental_var_name' at the specified timestep into the 'elem_var_value_map' which is passed in.

Definition at line 999 of file exodusII_io_helper.C.

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

1002 {
1003  this->read_var_names(ELEMENTAL);
1004 
1005  // See if we can find the variable we are looking for
1006  std::size_t var_index = 0;
1007  bool found = false;
1008 
1009  // Do a linear search for nodal_var_name in nodal_var_names
1010  for (; var_index<elem_var_names.size(); ++var_index)
1011  if (elem_var_names[var_index] == elemental_var_name)
1012  {
1013  found = true;
1014  break;
1015  }
1016 
1017  if (!found)
1018  {
1019  libMesh::err << "Available variables: " << std::endl;
1020  for (std::size_t i=0; i<elem_var_names.size(); ++i)
1021  libMesh::err << elem_var_names[i] << std::endl;
1022 
1023  libmesh_error_msg("Unable to locate variable named: " << elemental_var_name);
1024  }
1025 
1026  // Sequential index which we can use to look up the element ID in the elem_num_map.
1027  unsigned ex_el_num = 0;
1028 
1029  for (unsigned i=0; i<static_cast<unsigned>(num_elem_blk); i++)
1030  {
1032  block_ids[i],
1036  libmesh_nullptr);
1037  EX_CHECK_ERR(ex_err, "Error getting number of elements in block.");
1038 
1039  std::vector<Real> block_elem_var_values(num_elem);
1041  time_step,
1042  var_index+1,
1043  block_ids[i],
1045  &block_elem_var_values[0]);
1046  EX_CHECK_ERR(ex_err, "Error getting elemental values.");
1047 
1048  for (unsigned j=0; j<static_cast<unsigned>(num_elem_this_blk); j++)
1049  {
1050  // Use the elem_num_map to obtain the ID of this element in the Exodus file,
1051  // and remember to subtract 1 since libmesh is zero-based and Exodus is 1-based.
1052  unsigned mapped_elem_id = this->elem_num_map[ex_el_num] - 1;
1053 
1054  // Store the elemental value in the map.
1055  elem_var_value_map[mapped_elem_id] = block_elem_var_values[j];
1056 
1057  // Go to the next sequential element ID.
1058  ex_el_num++;
1059  }
1060  }
1061 }
EXODUS_EXPORT int ex_get_elem_block(int exoid, ex_entity_id elem_blk_id, char *elem_type, void_int *num_elem_this_blk, void_int *num_nodes_per_elem, void_int *num_attr)
std::vector< std::string > elem_var_names
const class libmesh_nullptr_t libmesh_nullptr
void read_var_names(ExodusVarType type)
OStreamProxy err(std::cerr)
EXODUS_EXPORT int ex_get_elem_var(int exoid, int time_step, int elem_var_index, ex_entity_id elem_blk_id, int64_t num_elem_this_blk, void *elem_var_vals)
void libMesh::ExodusII_IO_Helper::read_header ( )

Reads an ExodusII mesh file header.

Definition at line 365 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.

366 {
368  &title[0],
369  &num_dim,
370  &num_nodes,
371  &num_elem,
372  &num_elem_blk,
373  &num_node_sets,
374  &num_side_sets);
375 
376  EX_CHECK_ERR(ex_err, "Error retrieving header info.");
377 
378  this->read_num_time_steps();
379 
381  EX_CHECK_ERR(ex_err, "Error reading number of nodal variables.");
382 
384  EX_CHECK_ERR(ex_err, "Error reading number of elemental variables.");
385 
387  EX_CHECK_ERR(ex_err, "Error reading number of global variables.");
388 
389  message("Exodus header info retrieved successfully.");
390 }
EXODUS_EXPORT int ex_get_init(int exoid, char *title, void_int *num_dim, void_int *num_nodes, void_int *num_elem, void_int *num_elem_blk, void_int *num_node_sets, void_int *num_side_sets)
void message(const std::string &msg)
EXODUS_EXPORT int ex_get_var_param(int exoid, const char *var_type, int *num_vars)
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 843 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().

844 {
845  // Read the nodal variable names from file, so we can see if we have the one we're looking for
846  this->read_var_names(NODAL);
847 
848  // See if we can find the variable we are looking for
849  std::size_t var_index = 0;
850  bool found = false;
851 
852  // Do a linear search for nodal_var_name in nodal_var_names
853  for (; var_index<nodal_var_names.size(); ++var_index)
854  {
855  found = (nodal_var_names[var_index] == nodal_var_name);
856  if (found)
857  break;
858  }
859 
860  if (!found)
861  {
862  libMesh::err << "Available variables: " << std::endl;
863  for (std::size_t i=0; i<nodal_var_names.size(); ++i)
864  libMesh::err << nodal_var_names[i] << std::endl;
865 
866  libmesh_error_msg("Unable to locate variable named: " << nodal_var_name);
867  }
868 
869  // Allocate enough space to store the nodal variable values
870  nodal_var_values.resize(num_nodes);
871 
872  // Call the Exodus API to read the nodal variable values
874  time_step,
875  var_index+1,
876  num_nodes,
877  &nodal_var_values[0]);
878  EX_CHECK_ERR(ex_err, "Error reading nodal variable values!");
879 }
std::vector< Real > nodal_var_values
EXODUS_EXPORT int ex_get_nodal_var(int exoid, int time_step, int nodal_var_index, int64_t num_nodes, void *nodal_var_vals)
void read_var_names(ExodusVarType type)
OStreamProxy err(std::cerr)
std::vector< std::string > nodal_var_names
void libMesh::ExodusII_IO_Helper::read_node_num_map ( )

Reads the optional node_num_map from the ExodusII mesh file.

Definition at line 476 of file exodusII_io_helper.C.

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

477 {
478  node_num_map.resize(num_nodes);
479 
481  node_num_map.empty() ? libmesh_nullptr : &node_num_map[0]);
482 
483  EX_CHECK_ERR(ex_err, "Error retrieving nodal number map.");
484  message("Nodal numbering map retrieved successfully.");
485 
486  if (verbose)
487  {
488  libMesh::out << "[" << this->processor_id() << "] node_num_map[i] = ";
489  for (unsigned int i=0; i<static_cast<unsigned int>(std::min(10, num_nodes-1)); ++i)
490  libMesh::out << node_num_map[i] << ", ";
491  libMesh::out << "... " << node_num_map.back() << std::endl;
492  }
493 }
EXODUS_EXPORT int ex_get_node_num_map(int exoid, void_int *node_map)
const class libmesh_nullptr_t libmesh_nullptr
void message(const std::string &msg)
OStreamProxy out(std::cout)
long double min(long double a, double b)
processor_id_type processor_id() const
void libMesh::ExodusII_IO_Helper::read_nodes ( )

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

Definition at line 459 of file exodusII_io_helper.C.

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

460 {
461  x.resize(num_nodes);
462  y.resize(num_nodes);
463  z.resize(num_nodes);
464 
466  static_cast<void *>(&x[0]),
467  static_cast<void *>(&y[0]),
468  static_cast<void *>(&z[0]));
469 
470  EX_CHECK_ERR(ex_err, "Error retrieving nodal data.");
471  message("Nodal data retrieved successfully.");
472 }
EXODUS_EXPORT int ex_get_coord(int exoid, void *x_coor, void *y_coor, void *z_coor)
void message(const std::string &msg)
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 752 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.

753 {
754  libmesh_assert_less (static_cast<unsigned int>(id), nodeset_ids.size());
755  libmesh_assert_less (static_cast<unsigned int>(id), num_nodes_per_set.size());
756  libmesh_assert_less (static_cast<unsigned int>(id), num_node_df_per_set.size());
757 
759  nodeset_ids[id],
760  &num_nodes_per_set[id],
761  &num_node_df_per_set[id]);
762  EX_CHECK_ERR(ex_err, "Error retrieving nodeset parameters.");
763  message("Parameters retrieved successfully for nodeset: ", id);
764 
765  node_list.resize(num_nodes_per_set[id]);
766 
767  // Don't call ex_get_node_set unless there are actually nodes there to get.
768  // Exodus prints an annoying warning message in DEBUG mode otherwise...
769  if (num_nodes_per_set[id] > 0)
770  {
771  ex_err = exII::ex_get_node_set(ex_id,
772  nodeset_ids[id],
773  &node_list[0]);
774 
775  EX_CHECK_ERR(ex_err, "Error retrieving nodeset data.");
776  message("Data retrieved successfully for nodeset: ", id);
777  }
778 }
EXODUS_EXPORT int ex_get_node_set_param(int exoid, ex_entity_id node_set_id, void_int *num_nodes_in_set, void_int *num_df_in_set)
EXODUS_EXPORT int ex_get_node_set(int exoid, ex_entity_id node_set_id, void_int *node_set_node_list)
std::vector< int > num_node_df_per_set
void message(const std::string &msg)
std::vector< int > num_nodes_per_set
void libMesh::ExodusII_IO_Helper::read_nodeset_info ( )

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

Definition at line 681 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.

682 {
683  nodeset_ids.resize(num_node_sets);
684  if (num_node_sets > 0)
685  {
687  &nodeset_ids[0]);
688  EX_CHECK_ERR(ex_err, "Error retrieving nodeset information.");
689  message("All nodeset information retrieved successfully.");
690 
691  // Resize appropriate data structures -- only do this once outnode the loop
694  }
695 
696  char name_buffer[MAX_STR_LENGTH+1];
697  for (int i=0; i<num_node_sets; ++i)
698  {
700  nodeset_ids[i], name_buffer);
701  EX_CHECK_ERR(ex_err, "Error getting node set name.");
702  id_to_ns_names[nodeset_ids[i]] = name_buffer;
703  }
704  message("All node set names retrieved successfully.");
705 }
EXODUS_EXPORT int ex_get_name(int exoid, ex_entity_type obj_type, ex_entity_id entity_id, char *name)
EXODUS_EXPORT int ex_get_node_set_ids(int exoid, void_int *ids)
std::vector< int > num_node_df_per_set
void message(const std::string &msg)
std::map< int, std::string > id_to_ns_names
std::vector< int > num_nodes_per_set
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 835 of file exodusII_io_helper.C.

References EX_INQ_TIME, inquire(), and num_time_steps.

Referenced by read_header(), and read_time_steps().

836 {
838  this->inquire(exII::EX_INQ_TIME, "Error retrieving number of time steps");
839 }
int inquire(int req_info, std::string error_msg="")
void libMesh::ExodusII_IO_Helper::read_qa_records ( )

Reads the QA records from an ExodusII file. We can use this to detect when e.g. CUBIT 14+ was used to generate a Mesh file, and work around certain known bugs in that version.

Definition at line 395 of file exodusII_io_helper.C.

References ex_err, ex_get_qa(), ex_id, EX_INQ_QA, inquire(), libMesh::out, and verbose.

396 {
397  // The QA records are four MAX_STR_LENGTH-byte character strings.
398  int num_qa_rec =
399  this->inquire(exII::EX_INQ_QA, "Error retrieving number of QA records");
400 
401  if (verbose)
402  libMesh::out << "Found "
403  << num_qa_rec
404  << " QA record(s) in the Exodus file."
405  << std::endl;
406 
407  if (num_qa_rec > 0)
408  {
409  // How to dynamically allocate an array of fixed-size char * arrays in C++.
410  // http://stackoverflow.com/questions/8529359/creating-a-dynamic-sized-array-of-fixed-sized-int-arrays-in-c
411  typedef char * inner_array_t[4];
412  inner_array_t * qa_record = new inner_array_t[num_qa_rec];
413 
414  for (int i=0; i<num_qa_rec; i++)
415  for (int j=0; j<4; j++)
416  qa_record[i][j] = new char[MAX_STR_LENGTH+1];
417 
418  ex_err = exII::ex_get_qa (ex_id, qa_record);
419  EX_CHECK_ERR(ex_err, "Error reading the QA records.");
420 
421  // Print the QA records
422  if (verbose)
423  {
424  for (int i=0; i<num_qa_rec; i++)
425  {
426  libMesh::out << "QA Record: " << i << std::endl;
427  for (int j=0; j<4; j++)
428  libMesh::out << qa_record[i][j] << std::endl;
429  }
430  }
431 
432 
433  // Clean up dynamically-allocated memory
434  for (int i=0; i<num_qa_rec; i++)
435  for (int j=0; j<4; j++)
436  delete [] qa_record[i][j];
437 
438  delete [] qa_record;
439  }
440 }
EXODUS_EXPORT int ex_get_qa(int exoid, char *qa_record[][4])
int inquire(int req_info, std::string error_msg="")
OStreamProxy out(std::cout)
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 709 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.

710 {
711  libmesh_assert_less (static_cast<unsigned int>(id), ss_ids.size());
712  libmesh_assert_less (static_cast<unsigned int>(id), num_sides_per_set.size());
713  libmesh_assert_less (static_cast<unsigned int>(id), num_df_per_set.size());
714  libmesh_assert_less_equal (static_cast<unsigned int>(offset), elem_list.size());
715  libmesh_assert_less_equal (static_cast<unsigned int>(offset), side_list.size());
716 
718  ss_ids[id],
719  &num_sides_per_set[id],
720  &num_df_per_set[id]);
721  EX_CHECK_ERR(ex_err, "Error retrieving sideset parameters.");
722  message("Parameters retrieved successfully for sideset: ", id);
723 
724 
725  // It's OK for offset==elem_list.size() as long as num_sides_per_set[id]==0
726  // because in that case we don't actually read anything...
727 #ifdef DEBUG
728  if (static_cast<unsigned int>(offset) == elem_list.size() ||
729  static_cast<unsigned int>(offset) == side_list.size() )
730  libmesh_assert_equal_to (num_sides_per_set[id], 0);
731 #endif
732 
733 
734  // Don't call ex_get_side_set unless there are actually sides there to get.
735  // Exodus prints an annoying warning in DEBUG mode otherwise...
736  if (num_sides_per_set[id] > 0)
737  {
738  ex_err = exII::ex_get_side_set(ex_id,
739  ss_ids[id],
740  &elem_list[offset],
741  &side_list[offset]);
742  EX_CHECK_ERR(ex_err, "Error retrieving sideset data.");
743  message("Data retrieved successfully for sideset: ", id);
744 
745  for (int i=0; i<num_sides_per_set[id]; i++)
746  id_list[i+offset] = ss_ids[id];
747  }
748 }
std::vector< int > num_sides_per_set
EXODUS_EXPORT int ex_get_side_set(int exoid, ex_entity_id side_set_id, void_int *side_set_elem_list, void_int *side_set_side_list)
EXODUS_EXPORT int ex_get_side_set_param(int exoid, ex_entity_id side_set_id, void_int *num_side_in_set, void_int *num_dist_fact_in_set)
void message(const std::string &msg)
std::vector< int > num_df_per_set
void libMesh::ExodusII_IO_Helper::read_sideset_info ( )

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

Definition at line 648 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.

649 {
650  ss_ids.resize(num_side_sets);
651  if (num_side_sets > 0)
652  {
654  &ss_ids[0]);
655  EX_CHECK_ERR(ex_err, "Error retrieving sideset information.");
656  message("All sideset information retrieved successfully.");
657 
658  // Resize appropriate data structures -- only do this once outside the loop
661 
662  // Inquire about the length of the concatenated side sets element list
663  num_elem_all_sidesets = inquire(exII::EX_INQ_SS_ELEM_LEN, "Error retrieving length of the concatenated side sets element list!");
664 
668  }
669 
670  char name_buffer[MAX_STR_LENGTH+1];
671  for (int i=0; i<num_side_sets; ++i)
672  {
674  ss_ids[i], name_buffer);
675  EX_CHECK_ERR(ex_err, "Error getting side set name.");
676  id_to_ss_names[ss_ids[i]] = name_buffer;
677  }
678  message("All side set names retrieved successfully.");
679 }
EXODUS_EXPORT int ex_get_side_set_ids(int exoid, void_int *ids)
EXODUS_EXPORT int ex_get_name(int exoid, ex_entity_type obj_type, ex_entity_id entity_id, char *name)
std::vector< int > num_sides_per_set
std::map< int, std::string > id_to_ss_names
int inquire(int req_info, std::string error_msg="")
void message(const std::string &msg)
std::vector< int > num_df_per_set
void libMesh::ExodusII_IO_Helper::read_time_steps ( )

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

Definition at line 820 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.

821 {
822  // Make sure we have an up-to-date count of the number of time steps in the file.
823  this->read_num_time_steps();
824 
825  if (num_time_steps > 0)
826  {
827  time_steps.resize(num_time_steps);
829  EX_CHECK_ERR(ex_err, "Error reading timesteps!");
830  }
831 }
std::vector< Real > time_steps
EXODUS_EXPORT int ex_get_all_times(int exoid, void *time_values)
void libMesh::ExodusII_IO_Helper::read_var_names ( ExodusVarType  type)

Definition at line 883 of file exodusII_io_helper.C.

References elem_var_names, ELEMENTAL, 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().

884 {
885  switch (type)
886  {
887  case NODAL:
889  break;
890  case ELEMENTAL:
892  break;
893  case GLOBAL:
895  break;
896  default:
897  libmesh_error_msg("Unrecognized ExodusVarType " << type);
898  }
899 }
std::vector< std::string > elem_var_names
void read_var_names_impl(const char *var_type, int &count, std::vector< std::string > &result)
std::vector< std::string > global_var_names
std::vector< std::string > nodal_var_names
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 903 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().

906 {
907  // First read and store the number of names we have
908  ex_err = exII::ex_get_var_param(ex_id, var_type, &count);
909  EX_CHECK_ERR(ex_err, "Error reading number of variables.");
910 
911  // Do nothing if no variables are detected
912  if (count == 0)
913  return;
914 
915  // Second read the actual names and convert them into a format we can use
916  NamesData names_table(count, MAX_STR_LENGTH);
917 
919  var_type,
920  count,
921  names_table.get_char_star_star()
922  );
923  EX_CHECK_ERR(ex_err, "Error reading variable names!");
924 
925  if (verbose)
926  {
927  libMesh::out << "Read the variable(s) from the file:" << std::endl;
928  for (int i=0; i<count; i++)
929  libMesh::out << names_table.get_char_star(i) << std::endl;
930  }
931 
932  // Allocate enough space for our variable name strings.
933  result.resize(count);
934 
935  // Copy the char buffers into strings.
936  for (int i=0; i<count; i++)
937  result[i] = names_table.get_char_star(i); // calls string::op=(const char *)
938 }
EXODUS_EXPORT int ex_get_var_names(int exoid, const char *var_type, int num_vars, char *var_names[])
EXODUS_EXPORT int ex_get_var_param(int exoid, const char *var_type, int *num_vars)
OStreamProxy out(std::cout)
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 1976 of file exodusII_io_helper.C.

References _coordinate_offset.

1977 {
1978  _coordinate_offset = p;
1979 }
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 1962 of file exodusII_io_helper.C.

References _use_mesh_dimension_instead_of_spatial_dimension.

void libMesh::ExodusII_IO_Helper::write_as_dimension ( unsigned  dim)

Sets the value of _write_as_dimension.

This directly controls the num_dim which is written to the Exodus file. If non-zero, this value supersedes all other dimensions, including: 1.) MeshBase::spatial_dimension() 2.) MeshBase::mesh_dimension() 3.) Any value passed to use_mesh_dimension_instead_of_spatial_dimension() This is useful/necessary for working around a bug in Paraview which prevents the "Plot Over Line" filter from working on 1D meshes.

Definition at line 1969 of file exodusII_io_helper.C.

References _write_as_dimension.

1970 {
1971  _write_as_dimension = dim;
1972 }
void libMesh::ExodusII_IO_Helper::write_element_values ( const MeshBase mesh,
const std::vector< Real > &  values,
int  timestep 
)

Writes the vector of values to the element variables.

Definition at line 1804 of file exodusII_io_helper.C.

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

1805 {
1806  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1807  return;
1808 
1809  // Loop over the element blocks and write the data one block at a time
1810  std::map<unsigned int, std::vector<unsigned int> > subdomain_map;
1811 
1812  // Ask the file how many element vars it has, store it in the num_elem_vars variable.
1814  EX_CHECK_ERR(ex_err, "Error reading number of elemental variables.");
1815 
1818 
1819  // loop through element and map between block and element vector
1820  for (; mesh_it!=end; ++mesh_it)
1821  {
1822  const Elem * elem = *mesh_it;
1823  subdomain_map[elem->subdomain_id()].push_back(elem->id());
1824  }
1825 
1826  // Use mesh.n_elem() to access into the values vector rather than
1827  // the number of elements the Exodus writer thinks the mesh has,
1828  // which may not include inactive elements.
1829  dof_id_type n_elem = mesh.n_elem();
1830 
1831  // For each variable, create a 'data' array which holds all the elemental variable
1832  // values *for a given block* on this processor, then write that data vector to file
1833  // before moving onto the next block.
1834  for (unsigned int i=0; i<static_cast<unsigned>(num_elem_vars); ++i)
1835  {
1836  // The size of the subdomain map is the number of blocks.
1837  std::map<unsigned int, std::vector<unsigned int> >::iterator it = subdomain_map.begin();
1838 
1839  for (unsigned int j=0; it!=subdomain_map.end(); ++it, ++j)
1840  {
1841  const std::vector<unsigned int> & elem_nums = (*it).second;
1842  const unsigned int num_elems_this_block =
1843  cast_int<unsigned int>(elem_nums.size());
1844  std::vector<Real> data(num_elems_this_block);
1845 
1846  for (unsigned int k=0; k<num_elems_this_block; ++k)
1847  data[k] = values[i*n_elem + elem_nums[k]];
1848 
1849  if (_single_precision)
1850  {
1851  std::vector<float> cast_data(data.begin(), data.end());
1852 
1854  timestep,
1855  i+1,
1856  this->get_block_id(j),
1857  num_elems_this_block,
1858  &cast_data[0]);
1859  }
1860  else
1861  {
1862  ex_err = exII::ex_put_elem_var(ex_id,
1863  timestep,
1864  i+1,
1865  this->get_block_id(j),
1866  num_elems_this_block,
1867  &data[0]);
1868  }
1869  EX_CHECK_ERR(ex_err, "Error writing element values.");
1870  }
1871  }
1872 
1874  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
1875 }
dof_id_type n_elem(const MeshBase::const_element_iterator &begin, const MeshBase::const_element_iterator &end)
Definition: mesh_tools.C:631
The base class for all geometric element types.
Definition: elem.h:86
IterBase * end
EXODUS_EXPORT int ex_put_elem_var(int exoid, int time_step, int elem_var_index, ex_entity_id elem_blk_id, int64_t num_elem_this_blk, const void *elem_var_vals)
subdomain_id_type subdomain_id() const
Definition: elem.h:1733
virtual element_iterator active_elements_begin()=0
virtual element_iterator active_elements_end()=0
EXODUS_EXPORT int ex_get_var_param(int exoid, const char *var_type, int *num_vars)
EXODUS_EXPORT int ex_update(int exoid)
IterBase * data
dof_id_type id() const
Definition: dof_object.h:624
virtual dof_id_type n_elem() const =0
uint8_t dof_id_type
Definition: id_types.h:64
processor_id_type processor_id() const
void libMesh::ExodusII_IO_Helper::write_elements ( const MeshBase mesh,
bool  use_discontinuous = false 
)
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 1318 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, elem_num_map, libMesh::MeshBase::elem_ref(), end, libMesh::Utility::enum_to_string(), 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::infinite(), libmesh_elem_num_to_exodus, libmesh_node_num_to_exodus, libMesh::MeshBase::n_active_elem(), libMesh::Elem::n_nodes(), libMesh::Elem::node_id(), 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.

1319 {
1320  // n_active_elem() is a parallel_only function
1321  unsigned int n_active_elem = mesh.n_active_elem();
1322 
1323  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1324  return;
1325 
1326  // Map from block ID to a vector of element IDs in that block. Element
1327  // IDs are now of type dof_id_type, subdomain IDs are of type subdomain_id_type.
1328  typedef std::map<subdomain_id_type, std::vector<dof_id_type> > subdomain_map_type;
1329  subdomain_map_type subdomain_map;
1330 
1333  //loop through element and map between block and element vector
1334  for (; mesh_it!=end; ++mesh_it)
1335  {
1336  const Elem * elem = *mesh_it;
1337  subdomain_map[ elem->subdomain_id() ].push_back(elem->id());
1338  }
1339 
1340  // element map vector
1341  num_elem_blk = cast_int<int>(subdomain_map.size());
1342  block_ids.resize(num_elem_blk);
1343  elem_num_map.resize(n_active_elem);
1344  std::vector<int>::iterator curr_elem_map_end = elem_num_map.begin();
1345 
1346  // Note: It appears that there is a bug in exodusII::ex_put_name where
1347  // the index returned from the ex_id_lkup is erronously used. For now
1348  // the work around is to use the alternative function ex_put_names, but
1349  // this function requires a char ** datastructure.
1350  NamesData names_table(num_elem_blk, MAX_STR_LENGTH);
1351 
1352  // This counter is used to fill up the libmesh_elem_num_to_exodus map in the loop below.
1353  unsigned libmesh_elem_num_to_exodus_counter = 0;
1354 
1355  // counter indexes into the block_ids vector
1356  unsigned int counter = 0;
1357 
1358  // node counter for discontinuous plotting
1359  unsigned int node_counter = 0;
1360 
1361  for (subdomain_map_type::iterator it=subdomain_map.begin(); it!=subdomain_map.end(); ++it)
1362  {
1363  block_ids[counter] = (*it).first;
1364  names_table.push_back_entry(mesh.subdomain_name((*it).first));
1365 
1366  // Get a reference to a vector of element IDs for this subdomain.
1367  subdomain_map_type::mapped_type & tmp_vec = (*it).second;
1368 
1370 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1371  // Skip infinite element-blocks; they can not be viewed in most visualization software
1372  // as paraview.
1373  if (mesh.elem_ref(tmp_vec[0]).infinite())
1374  continue;
1375 #endif
1376 
1377  //Use the first element in this block to get representative information.
1378  //Note that Exodus assumes all elements in a block are of the same type!
1379  //We are using that same assumption here!
1380  const ExodusII_IO_Helper::Conversion conv =
1381  em.assign_conversion(mesh.elem_ref(tmp_vec[0]).type());
1382  num_nodes_per_elem = mesh.elem_ref(tmp_vec[0]).n_nodes();
1383 
1385  (*it).first,
1386  conv.exodus_elem_type().c_str(),
1387  tmp_vec.size(),
1389  /*num_attr=*/0);
1390 
1391  EX_CHECK_ERR(ex_err, "Error writing element block.");
1392 
1393  connect.resize(tmp_vec.size()*num_nodes_per_elem);
1394 
1395  for (std::size_t i=0; i<tmp_vec.size(); i++)
1396  {
1397  unsigned int elem_id = tmp_vec[i];
1398  libmesh_elem_num_to_exodus[elem_id] = ++libmesh_elem_num_to_exodus_counter; // 1-based indexing for Exodus
1399 
1400  const Elem & elem = mesh.elem_ref(elem_id);
1401 
1402  // We *might* be able to get away with writing mixed element
1403  // types which happen to have the same number of nodes, but
1404  // do we actually *want* to get away with that?
1405  // .) No visualization software would be able to handle it.
1406  // .) There'd be no way for us to read it back in reliably.
1407  // .) Even elements with the same number of nodes may have different connectivities (?)
1408 
1409  // This needs to be more than an assert so we don't fail
1410  // with a mysterious segfault while trying to write mixed
1411  // element meshes in optimized mode.
1412  if (elem.type() != conv.get_canonical_type())
1413  libmesh_error_msg("Error: Exodus requires all elements with a given subdomain ID to be the same type.\n" \
1414  << "Can't write both " \
1415  << Utility::enum_to_string(elem.type()) \
1416  << " and " \
1418  << " in the same block!");
1419 
1420 
1421  for (unsigned int j=0; j<static_cast<unsigned int>(num_nodes_per_elem); ++j)
1422  {
1423  unsigned connect_index = (i*num_nodes_per_elem)+j;
1424  unsigned elem_node_index = conv.get_inverse_node_map(j); // inverse node map is for writing.
1425  if (verbose)
1426  {
1427  libMesh::out << "Exodus node index " << j
1428  << " = LibMesh node index " << elem_node_index << std::endl;
1429  }
1430 
1431  if (!use_discontinuous)
1432  {
1433  // The global id for the current node in libmesh.
1434  dof_id_type libmesh_node_id = elem.node_id(elem_node_index);
1435 
1436  // Find the zero-based libmesh id in the map, this
1437  // should be faster than doing linear searches on
1438  // the node_num_map.
1439  std::map<int, int>::iterator pos =
1440  libmesh_node_num_to_exodus.find(cast_int<int>(libmesh_node_id));
1441 
1442  // Make sure it was found.
1443  if (pos == libmesh_node_num_to_exodus.end())
1444  libmesh_error_msg("libmesh node id " << libmesh_node_id << " not found in node_num_map.");
1445 
1446  // Write the Exodus global node id associated with
1447  // this libmesh node number to the connectivity
1448  // array.
1449  connect[connect_index] = pos->second;
1450  }
1451  else
1452  {
1453  // FIXME: We are hard-coding the 1-based node
1454  // numbering assumption here, so writing
1455  // "discontinuous" Exodus files won't work with node
1456  // numberings that have "holes".
1457  connect[connect_index] = node_counter + elem_node_index + 1;
1458  }
1459  }
1460 
1461  node_counter += num_nodes_per_elem;
1462  }
1463 
1464  ex_err = exII::ex_put_elem_conn(ex_id, (*it).first, &connect[0]);
1465  EX_CHECK_ERR(ex_err, "Error writing element connectivities");
1466 
1467  // This transform command stores its result in a range that begins at the third argument,
1468  // so this command is adding values to the elem_num_map vector starting from curr_elem_map_end.
1469  curr_elem_map_end = std::transform(tmp_vec.begin(),
1470  tmp_vec.end(),
1471  curr_elem_map_end,
1472  std::bind2nd(std::plus<subdomain_map_type::mapped_type::value_type>(), 1)); // Adds one to each id to make a 1-based exodus file!
1473 
1474  // But if we don't want to add one, we just want to put the values
1475  // of tmp_vec into elem_map in the right location, we can use
1476  // std::copy().
1477  // curr_elem_map_end = std::copy(tmp_vec.begin(), tmp_vec.end(), curr_elem_map_end);
1478 
1479  counter++;
1480  }
1481 
1482  // write out the element number map that we created
1484  EX_CHECK_ERR(ex_err, "Error writing element map");
1485 
1486  // Write out the block names
1487  if (num_elem_blk > 0)
1488  {
1489  ex_err = exII::ex_put_names(ex_id, exII::EX_ELEM_BLOCK, names_table.get_char_star_star());
1490  EX_CHECK_ERR(ex_err, "Error writing element names");
1491  }
1492 }
virtual dof_id_type n_active_elem() const =0
virtual ElemType type() const =0
The base class for all geometric element types.
Definition: elem.h:86
ExodusII_IO_Helper::Conversion assign_conversion(std::string type_str)
EXODUS_EXPORT int ex_put_elem_block(int exoid, ex_entity_id elem_blk_id, const char *elem_type, int64_t num_elem_this_blk, int64_t num_nodes_per_elem, int64_t num_attr)
IterBase * end
std::map< int, int > libmesh_elem_num_to_exodus
virtual bool infinite() const =0
virtual unsigned int n_nodes() const =0
std::string & subdomain_name(subdomain_id_type id)
Definition: mesh_base.C:573
subdomain_id_type subdomain_id() const
Definition: elem.h:1733
EXODUS_EXPORT int ex_put_names(int exoid, ex_entity_type obj_type, char *names[])
std::map< int, int > libmesh_node_num_to_exodus
virtual element_iterator active_elements_begin()=0
std::string enum_to_string(const T e)
virtual element_iterator active_elements_end()=0
virtual const Elem & elem_ref(const dof_id_type i) const
Definition: mesh_base.h:453
EXODUS_EXPORT int ex_put_elem_conn(int exoid, ex_entity_id elem_blk_id, const void_int *connect)
dof_id_type node_id(const unsigned int i) const
Definition: elem.h:1617
dof_id_type id() const
Definition: dof_object.h:624
EXODUS_EXPORT int ex_put_elem_num_map(int exoid, const void_int *elem_map)
OStreamProxy out(std::cout)
uint8_t dof_id_type
Definition: id_types.h:64
processor_id_type processor_id() const
void libMesh::ExodusII_IO_Helper::write_global_values ( const std::vector< Real > &  values,
int  timestep 
)

Writes the vector of global variables.

Definition at line 1940 of file exodusII_io_helper.C.

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

1941 {
1942  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1943  return;
1944 
1945  if (_single_precision)
1946  {
1947  std::vector<float> cast_values(values.begin(), values.end());
1948  ex_err = exII::ex_put_glob_vars(ex_id, timestep, num_global_vars, &cast_values[0]);
1949  }
1950  else
1951  {
1952  ex_err = exII::ex_put_glob_vars(ex_id, timestep, num_global_vars, &values[0]);
1953  }
1954  EX_CHECK_ERR(ex_err, "Error writing global values.");
1955 
1957  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
1958 }
EXODUS_EXPORT int ex_put_glob_vars(int exoid, int time_step, int num_glob_vars, const void *glob_var_vals)
EXODUS_EXPORT int ex_update(int exoid)
processor_id_type processor_id() const
void libMesh::ExodusII_IO_Helper::write_information_records ( const std::vector< std::string > &  records)

Writes the vector of information records.

Definition at line 1901 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().

1902 {
1903  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1904  return;
1905 
1906  // There may already be information records in the file (for
1907  // example, if we're appending) and in that case, according to the
1908  // Exodus documentation, writing more information records is not
1909  // supported.
1910  int num_info = inquire(exII::EX_INQ_INFO, "Error retrieving the number of information records from file!");
1911  if (num_info > 0)
1912  {
1913  libMesh::err << "Warning! The Exodus file already contains information records.\n"
1914  << "Exodus does not support writing additional records in this situation."
1915  << std::endl;
1916  return;
1917  }
1918 
1919  int num_records = cast_int<int>(records.size());
1920 
1921  if (num_records > 0)
1922  {
1923  NamesData info(num_records, MAX_LINE_LENGTH);
1924 
1925  // If an entry is longer than MAX_LINE_LENGTH characters it's not an error, we just
1926  // write the first MAX_LINE_LENGTH characters to the file.
1927  for (std::size_t i=0; i<records.size(); ++i)
1928  info.push_back_entry(records[i]);
1929 
1930  ex_err = exII::ex_put_info(ex_id, num_records, info.get_char_star_star());
1931  EX_CHECK_ERR(ex_err, "Error writing global values.");
1932 
1934  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
1935  }
1936 }
EXODUS_EXPORT int ex_put_info(int exoid, int num_info, char *info[])
int inquire(int req_info, std::string error_msg="")
OStreamProxy err(std::cerr)
EXODUS_EXPORT int ex_update(int exoid)
processor_id_type processor_id() const
void libMesh::ExodusII_IO_Helper::write_nodal_coordinates ( const MeshBase mesh,
bool  use_discontinuous = false 
)
virtual

Writes the nodal coordinates contained in "mesh"

Reimplemented in libMesh::Nemesis_IO_Helper.

Definition at line 1193 of file exodusII_io_helper.C.

References _coordinate_offset, _run_only_on_proc0, _single_precision, libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), end, ex_err, ex_id, ex_put_coord(), ex_put_node_num_map(), libMesh::DofObject::id(), libmesh_node_num_to_exodus, libmesh_nullptr, libMesh::Elem::n_nodes(), node_num_map, libMesh::MeshBase::nodes_begin(), libMesh::MeshBase::nodes_end(), num_nodes, libMesh::Elem::point(), libMesh::ParallelObject::processor_id(), x, y, and z.

1194 {
1195  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1196  return;
1197 
1198  // Clear existing data from any previous calls.
1199  x.clear();
1200  y.clear();
1201  z.clear();
1202  node_num_map.clear();
1203 
1204  // Reserve space in the nodal coordinate vectors. num_nodes is
1205  // exact, this just allows us to do away with one potentially
1206  // error-inducing loop index.
1207  x.reserve(num_nodes);
1208  y.reserve(num_nodes);
1209  z.reserve(num_nodes);
1210 
1211  // And in the node_num_map - since the nodes aren't organized in
1212  // blocks, libmesh will always write out the identity map
1213  // here... unless there has been some refinement and coarsening, or
1214  // node deletion without a corresponding call to contract(). You
1215  // need to write this any time there could be 'holes' in the node
1216  // numbering, so we write it every time.
1217  node_num_map.reserve(num_nodes);
1218 
1219  // Clear out any previously-mapped node IDs.
1221 
1222  if (!use_discontinuous)
1223  {
1226  for (; it != end; ++it)
1227  {
1228  const Node & node = **it;
1229 
1230  x.push_back(node(0) + _coordinate_offset(0));
1231 
1232 #if LIBMESH_DIM > 1
1233  y.push_back(node(1) + _coordinate_offset(1));
1234 #else
1235  y.push_back(0.);
1236 #endif
1237 #if LIBMESH_DIM > 2
1238  z.push_back(node(2) + _coordinate_offset(2));
1239 #else
1240  z.push_back(0.);
1241 #endif
1242 
1243  // Fill in node_num_map entry with the proper (1-based) node id
1244  node_num_map.push_back(node.id() + 1);
1245 
1246  // Also map the zero-based libmesh node id to the 1-based
1247  // Exodus ID it will be assigned (this is equivalent to the
1248  // current size of the x vector).
1249  libmesh_node_num_to_exodus[ cast_int<int>(node.id()) ] = cast_int<int>(x.size());
1250  }
1251  }
1252  else
1253  {
1256 
1257  for (; it!=end; ++it)
1258  {
1259  const Elem * elem = *it;
1260 
1261  for (unsigned int n=0; n<elem->n_nodes(); n++)
1262  {
1263  x.push_back(elem->point(n)(0));
1264 #if LIBMESH_DIM > 1
1265  y.push_back(elem->point(n)(1));
1266 #else
1267  y.push_back(0.);
1268 #endif
1269 #if LIBMESH_DIM > 2
1270  z.push_back(elem->point(n)(2));
1271 #else
1272  z.push_back(0.);
1273 #endif
1274 
1275  // Let's skip the node_num_map in the discontinuous
1276  // case, since we're effectively duplicating nodes for
1277  // the sake of discontinuous visualization, so it isn't
1278  // clear how to deal with node_num_map here. This means
1279  // that writing discontinuous meshes won't work with
1280  // element numberings that have "holes".
1281  }
1282  }
1283  }
1284 
1285  if (_single_precision)
1286  {
1287  std::vector<float>
1288  x_single(x.begin(), x.end()),
1289  y_single(y.begin(), y.end()),
1290  z_single(z.begin(), z.end());
1291 
1293  x_single.empty() ? libmesh_nullptr : &x_single[0],
1294  y_single.empty() ? libmesh_nullptr : &y_single[0],
1295  z_single.empty() ? libmesh_nullptr : &z_single[0]);
1296  }
1297  else
1298  {
1299  ex_err = exII::ex_put_coord(ex_id,
1300  x.empty() ? libmesh_nullptr : &x[0],
1301  y.empty() ? libmesh_nullptr : &y[0],
1302  z.empty() ? libmesh_nullptr : &z[0]);
1303  }
1304 
1305 
1306  EX_CHECK_ERR(ex_err, "Error writing coordinates to Exodus file.");
1307 
1308  if (!use_discontinuous)
1309  {
1310  // Also write the (1-based) node_num_map to the file.
1312  EX_CHECK_ERR(ex_err, "Error writing node_num_map");
1313  }
1314 }
A geometric point in (x,y,z) space associated with a DOF.
Definition: node.h:52
The base class for all geometric element types.
Definition: elem.h:86
const class libmesh_nullptr_t libmesh_nullptr
EXODUS_EXPORT int ex_put_coord(int exoid, const void *x_coor, const void *y_coor, const void *z_coor)
IterBase * end
virtual unsigned int n_nodes() const =0
virtual node_iterator nodes_begin()=0
std::map< int, int > libmesh_node_num_to_exodus
virtual element_iterator active_elements_begin()=0
virtual element_iterator active_elements_end()=0
virtual node_iterator nodes_end()=0
const Point & point(const unsigned int i) const
Definition: elem.h:1595
EXODUS_EXPORT int ex_put_node_num_map(int exoid, const void_int *node_map)
dof_id_type id() const
Definition: dof_object.h:624
processor_id_type processor_id() const
void libMesh::ExodusII_IO_Helper::write_nodal_values ( int  var_id,
const std::vector< Real > &  values,
int  timestep 
)

Writes the vector of values to a nodal variable.

Definition at line 1879 of file exodusII_io_helper.C.

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

Referenced by libMesh::Nemesis_IO_Helper::write_nodal_solution().

1880 {
1881  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1882  return;
1883 
1884  if (_single_precision)
1885  {
1886  std::vector<float> cast_values(values.begin(), values.end());
1887  ex_err = exII::ex_put_nodal_var(ex_id, timestep, var_id, num_nodes, &cast_values[0]);
1888  }
1889  else
1890  {
1891  ex_err = exII::ex_put_nodal_var(ex_id, timestep, var_id, num_nodes, &values[0]);
1892  }
1893  EX_CHECK_ERR(ex_err, "Error writing nodal values.");
1894 
1896  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
1897 }
EXODUS_EXPORT int ex_update(int exoid)
processor_id_type processor_id() const
EXODUS_EXPORT int ex_put_nodal_var(int exoid, int time_step, int nodal_var_index, int64_t num_nodes, const void *nodal_var_vals)
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 1612 of file exodusII_io_helper.C.

References _run_only_on_proc0, libMesh::BoundaryInfo::build_node_boundary_ids(), libMesh::BoundaryInfo::build_node_list(), ex_err, ex_id, EX_NODE_SET, ex_put_names(), ex_put_node_set(), ex_put_node_set_param(), libMesh::MeshBase::get_boundary_info(), libMesh::BoundaryInfo::get_nodeset_name(), libMesh::ParallelObject::processor_id(), and libMesh::ExodusII_IO_Helper::NamesData::push_back_entry().

1613 {
1614  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1615  return;
1616 
1617  std::vector< dof_id_type > nl;
1618  std::vector< boundary_id_type > il;
1619 
1620  mesh.get_boundary_info().build_node_list(nl, il);
1621 
1622  // Maps from nodeset id to the nodes
1623  std::map<boundary_id_type, std::vector<int> > node;
1624 
1625  // Accumulate the vectors to pass into ex_put_node_set
1626  for (std::size_t i=0; i<nl.size(); i++)
1627  node[il[i]].push_back(nl[i]+1);
1628 
1629  std::vector<boundary_id_type> node_boundary_ids;
1630  mesh.get_boundary_info().build_node_boundary_ids(node_boundary_ids);
1631 
1632  // Write out the nodeset names, but only if there is something to write
1633  if (node_boundary_ids.size() > 0)
1634  {
1635  NamesData names_table(node_boundary_ids.size(), MAX_STR_LENGTH);
1636 
1637  for (std::size_t i=0; i<node_boundary_ids.size(); i++)
1638  {
1639  boundary_id_type nodeset_id = node_boundary_ids[i];
1640 
1641  int actual_id = nodeset_id;
1642 
1643  names_table.push_back_entry(mesh.get_boundary_info().get_nodeset_name(nodeset_id));
1644 
1645  ex_err = exII::ex_put_node_set_param(ex_id, actual_id, node[nodeset_id].size(), 0);
1646  EX_CHECK_ERR(ex_err, "Error writing nodeset parameters");
1647 
1648  ex_err = exII::ex_put_node_set(ex_id, actual_id, &node[nodeset_id][0]);
1649  EX_CHECK_ERR(ex_err, "Error writing nodesets");
1650  }
1651 
1652  // Write out the nodeset names
1653  ex_err = exII::ex_put_names(ex_id, exII::EX_NODE_SET, names_table.get_char_star_star());
1654  EX_CHECK_ERR(ex_err, "Error writing nodeset names");
1655  }
1656 }
const BoundaryInfo & get_boundary_info() const
Definition: mesh_base.h:111
EXODUS_EXPORT int ex_put_node_set_param(int exoid, ex_entity_id node_set_id, int64_t num_nodes_in_set, int64_t num_dist_in_set)
const std::string & get_nodeset_name(boundary_id_type id) const
EXODUS_EXPORT int ex_put_node_set(int exoid, ex_entity_id node_set_id, const void_int *node_set_node_list)
int8_t boundary_id_type
Definition: id_types.h:51
EXODUS_EXPORT int ex_put_names(int exoid, ex_entity_type obj_type, char *names[])
void build_node_list(std::vector< dof_id_type > &node_id_list, std::vector< boundary_id_type > &bc_id_list) const
void build_node_boundary_ids(std::vector< boundary_id_type > &b_ids) const
processor_id_type processor_id() const
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

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 1497 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::BoundaryInfo::build_shellface_boundary_ids(), libMesh::BoundaryInfo::build_shellface_list(), libMesh::BoundaryInfo::build_side_boundary_ids(), libMesh::BoundaryInfo::build_side_list(), libMesh::MeshBase::elem_ptr(), libMesh::MeshBase::elem_ref(), ex_err, ex_id, ex_put_names(), ex_put_side_set(), ex_put_side_set_param(), EX_SIDE_SET, libMesh::MeshBase::get_boundary_info(), libMesh::ExodusII_IO_Helper::Conversion::get_inverse_shellface_map(), libMesh::ExodusII_IO_Helper::Conversion::get_inverse_side_map(), libMesh::BoundaryInfo::get_sideset_name(), libmesh_elem_num_to_exodus, libMesh::ParallelObject::processor_id(), libMesh::ExodusII_IO_Helper::NamesData::push_back_entry(), and side.

1498 {
1499  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1500  return;
1501 
1503 
1504  // Maps from sideset id to the element and sides
1505  std::map<int, std::vector<int> > elem;
1506  std::map<int, std::vector<int> > side;
1507  std::vector<boundary_id_type> side_boundary_ids;
1508 
1509  {
1510  std::vector< dof_id_type > el;
1511  std::vector< unsigned short int > sl;
1512  std::vector< boundary_id_type > il;
1513 
1514  mesh.get_boundary_info().build_side_list(el, sl, il);
1515 
1516  // Accumulate the vectors to pass into ex_put_side_set
1517  for (std::size_t i=0; i<el.size(); i++)
1518  {
1519  std::vector<const Elem *> family;
1520 #ifdef LIBMESH_ENABLE_AMR
1521 
1525  mesh.elem_ref(el[i]).active_family_tree_by_side(family, sl[i], false);
1526 #else
1527  family.push_back(mesh.elem_ptr(el[i]));
1528 #endif
1529 
1530  for (std::size_t j=0; j<family.size(); ++j)
1531  {
1532  const ExodusII_IO_Helper::Conversion conv =
1533  em.assign_conversion(mesh.elem_ptr(family[j]->id())->type());
1534 
1535  // Use the libmesh to exodus datastructure map to get the proper sideset IDs
1536  // The datastructure contains the "collapsed" contiguous ids
1537  elem[il[i]].push_back(libmesh_elem_num_to_exodus[family[j]->id()]);
1538  side[il[i]].push_back(conv.get_inverse_side_map(sl[i]));
1539  }
1540  }
1541 
1542  mesh.get_boundary_info().build_side_boundary_ids(side_boundary_ids);
1543  }
1544 
1545  {
1546  // add data for shell faces, if needed
1547 
1548  std::vector< dof_id_type > el;
1549  std::vector< unsigned short int > sl;
1550  std::vector< boundary_id_type > il;
1551 
1552  mesh.get_boundary_info().build_shellface_list(el, sl, il);
1553 
1554  // Accumulate the vectors to pass into ex_put_side_set
1555  for (std::size_t i=0; i<el.size(); i++)
1556  {
1557  std::vector<const Elem *> family;
1558 #ifdef LIBMESH_ENABLE_AMR
1559 
1563  mesh.elem_ref(el[i]).active_family_tree_by_side(family, sl[i], false);
1564 #else
1565  family.push_back(mesh.elem_ptr(el[i]));
1566 #endif
1567 
1568  for (std::size_t j=0; j<family.size(); ++j)
1569  {
1570  const ExodusII_IO_Helper::Conversion conv =
1571  em.assign_conversion(mesh.elem_ptr(family[j]->id())->type());
1572 
1573  // Use the libmesh to exodus datastructure map to get the proper sideset IDs
1574  // The datastructure contains the "collapsed" contiguous ids
1575  elem[il[i]].push_back(libmesh_elem_num_to_exodus[family[j]->id()]);
1576  side[il[i]].push_back(conv.get_inverse_shellface_map(sl[i]));
1577  }
1578  }
1579 
1580  std::vector<boundary_id_type> shellface_boundary_ids;
1581  mesh.get_boundary_info().build_shellface_boundary_ids(shellface_boundary_ids);
1582  for (std::size_t i=0; i<shellface_boundary_ids.size(); i++)
1583  side_boundary_ids.push_back(shellface_boundary_ids[i]);
1584  }
1585 
1586  // Write out the sideset names, but only if there is something to write
1587  if (side_boundary_ids.size() > 0)
1588  {
1589  NamesData names_table(side_boundary_ids.size(), MAX_STR_LENGTH);
1590 
1591  for (std::size_t i=0; i<side_boundary_ids.size(); i++)
1592  {
1593  boundary_id_type ss_id = side_boundary_ids[i];
1594  int actual_id = ss_id;
1595 
1596  names_table.push_back_entry(mesh.get_boundary_info().get_sideset_name(ss_id));
1597 
1598  ex_err = exII::ex_put_side_set_param(ex_id, actual_id, elem[ss_id].size(), 0);
1599  EX_CHECK_ERR(ex_err, "Error writing sideset parameters");
1600 
1601  ex_err = exII::ex_put_side_set(ex_id, actual_id, &elem[ss_id][0], &side[ss_id][0]);
1602  EX_CHECK_ERR(ex_err, "Error writing sidesets");
1603  }
1604 
1605  ex_err = exII::ex_put_names(ex_id, exII::EX_SIDE_SET, names_table.get_char_star_star());
1606  EX_CHECK_ERR(ex_err, "Error writing sideset names");
1607  }
1608 }
const BoundaryInfo & get_boundary_info() const
Definition: mesh_base.h:111
unsigned short int side
Definition: xdr_io.C:49
ExodusII_IO_Helper::Conversion assign_conversion(std::string type_str)
void build_side_list(std::vector< dof_id_type > &element_id_list, std::vector< unsigned short int > &side_list, std::vector< boundary_id_type > &bc_id_list) const
std::map< int, int > libmesh_elem_num_to_exodus
int8_t boundary_id_type
Definition: id_types.h:51
void build_side_boundary_ids(std::vector< boundary_id_type > &b_ids) const
EXODUS_EXPORT int ex_put_side_set(int exoid, ex_entity_id side_set_id, const void_int *side_set_elem_list, const void_int *side_set_side_list)
EXODUS_EXPORT int ex_put_names(int exoid, ex_entity_type obj_type, char *names[])
virtual const Elem & elem_ref(const dof_id_type i) const
Definition: mesh_base.h:453
EXODUS_EXPORT int ex_put_side_set_param(int exoid, ex_entity_id side_set_id, int64_t num_side_in_set, int64_t num_dist_fact_in_set)
void active_family_tree_by_side(std::vector< const Elem * > &family, const unsigned int side, const bool reset=true) const
Definition: elem.C:1671
void build_shellface_boundary_ids(std::vector< boundary_id_type > &b_ids) const
virtual const Elem * elem_ptr(const dof_id_type i) const =0
const std::string & get_sideset_name(boundary_id_type id) const
void build_shellface_list(std::vector< dof_id_type > &element_id_list, std::vector< unsigned short int > &shellface_list, std::vector< boundary_id_type > &bc_id_list) const
processor_id_type processor_id() const
void libMesh::ExodusII_IO_Helper::write_timestep ( int  timestep,
Real  time 
)

Writes the time for the timestep

Definition at line 1790 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().

1791 {
1792  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1793  return;
1794 
1795  ex_err = exII::ex_put_time(ex_id, timestep, &time);
1796  EX_CHECK_ERR(ex_err, "Error writing timestep.");
1797 
1799  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
1800 }
EXODUS_EXPORT int ex_update(int exoid)
EXODUS_EXPORT int ex_put_time(int exoid, int time_step, const void *time_value)
processor_id_type processor_id() const
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 943 of file exodusII_io_helper.C.

References ELEMENTAL, 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().

944 {
945  switch (type)
946  {
947  case NODAL:
948  this->write_var_names_impl("n", num_nodal_vars, names);
949  break;
950  case ELEMENTAL:
951  this->write_var_names_impl("e", num_elem_vars, names);
952  break;
953  case GLOBAL:
954  this->write_var_names_impl("g", num_global_vars, names);
955  break;
956  default:
957  libmesh_error_msg("Unrecognized ExodusVarType " << type);
958  }
959 }
void write_var_names_impl(const char *var_type, int &count, std::vector< std::string > &names)
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 963 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().

964 {
965  // Update the count variable so that it's available to other parts of the class.
966  count = cast_int<int>(names.size());
967 
968  // Write that number of variables to the file.
969  ex_err = exII::ex_put_var_param(ex_id, var_type, count);
970  EX_CHECK_ERR(ex_err, "Error setting number of vars.");
971 
972  if (names.size() > 0)
973  {
974  NamesData names_table(names.size(), MAX_STR_LENGTH);
975 
976  // Store the input names in the format required by Exodus.
977  for (std::size_t i=0; i<names.size(); ++i)
978  names_table.push_back_entry(names[i]);
979 
980  if (verbose)
981  {
982  libMesh::out << "Writing variable name(s) to file: " << std::endl;
983  for (std::size_t i=0; i<names.size(); ++i)
984  libMesh::out << names_table.get_char_star(i) << std::endl;
985  }
986 
987  ex_err = exII::ex_put_var_names(ex_id,
988  var_type,
989  cast_int<int>(names.size()),
990  names_table.get_char_star_star()
991  );
992 
993  EX_CHECK_ERR(ex_err, "Error writing variable names.");
994  }
995 }
EXODUS_EXPORT int ex_put_var_param(int exoid, const char *var_type, int num_vars)
EXODUS_EXPORT int ex_put_var_names(int exoid, const char *var_type, int num_vars, char *var_names[])
OStreamProxy out(std::cout)

Member Data Documentation

Point libMesh::ExodusII_IO_Helper::_coordinate_offset
protected

Definition at line 596 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 578 of file exodusII_io_helper.h.

Referenced by initialize_element_variables().

bool libMesh::ExodusII_IO_Helper::_global_vars_initialized
protected

Definition at line 581 of file exodusII_io_helper.h.

Referenced by initialize_global_variables().

bool libMesh::ExodusII_IO_Helper::_nodal_vars_initialized
protected

Definition at line 584 of file exodusII_io_helper.h.

Referenced by initialize_nodal_variables().

bool libMesh::ExodusII_IO_Helper::_use_mesh_dimension_instead_of_spatial_dimension
protected
unsigned libMesh::ExodusII_IO_Helper::_write_as_dimension
protected

Definition at line 593 of file exodusII_io_helper.h.

Referenced by initialize(), and write_as_dimension().

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

Definition at line 446 of file exodusII_io_helper.h.

Referenced by read_elem_in_block(), and write_elements().

std::string libMesh::ExodusII_IO_Helper::current_filename

Definition at line 560 of file exodusII_io_helper.h.

Referenced by create(), and open().

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

Definition at line 467 of file exodusII_io_helper.h.

Referenced by read_sideset(), and read_sideset_info().

std::vector<int> libMesh::ExodusII_IO_Helper::elem_num_map
std::vector<char> libMesh::ExodusII_IO_Helper::elem_type

Definition at line 497 of file exodusII_io_helper.h.

Referenced by ExodusII_IO_Helper(), get_elem_type(), and read_elem_in_block().

std::vector<std::string> libMesh::ExodusII_IO_Helper::elem_var_names
std::vector<Real> libMesh::ExodusII_IO_Helper::elem_var_values

Definition at line 531 of file exodusII_io_helper.h.

int libMesh::ExodusII_IO_Helper::ex_id

Definition at line 404 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_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_qa_records(), read_sideset(), read_sideset_info(), read_time_steps(), read_var_names_impl(), write_element_values(), libMesh::Nemesis_IO_Helper::write_elements(), write_elements(), libMesh::Nemesis_IO_Helper::write_exodus_initialization_info(), write_global_values(), write_information_records(), libMesh::Nemesis_IO_Helper::write_nodal_coordinates(), write_nodal_coordinates(), write_nodal_values(), libMesh::Nemesis_IO_Helper::write_nodesets(), write_nodesets(), libMesh::Nemesis_IO_Helper::write_sidesets(), 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 534 of file exodusII_io_helper.h.

Referenced by initialize_global_variables(), and read_var_names().

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

Definition at line 476 of file exodusII_io_helper.h.

Referenced by read_sideset(), and read_sideset_info().

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

Definition at line 537 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 539 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 538 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

Definition at line 522 of file exodusII_io_helper.h.

Referenced by read_nodal_var_values().

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

Definition at line 473 of file exodusII_io_helper.h.

Referenced by read_nodeset().

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

Definition at line 479 of file exodusII_io_helper.h.

Referenced by read_node_num_map(), and write_nodal_coordinates().

std::vector<int> libMesh::ExodusII_IO_Helper::nodeset_ids
int libMesh::ExodusII_IO_Helper::num_attr

Definition at line 437 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 461 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 440 of file exodusII_io_helper.h.

Referenced by read_sideset_info().

int libMesh::ExodusII_IO_Helper::num_elem_this_blk

Definition at line 431 of file exodusII_io_helper.h.

Referenced by read_elem_in_block(), and read_elemental_var_values().

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 464 of file exodusII_io_helper.h.

Referenced by read_nodeset(), and read_nodeset_info().

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

Definition at line 458 of file exodusII_io_helper.h.

Referenced by read_nodeset(), and read_nodeset_info().

int libMesh::ExodusII_IO_Helper::num_side_sets
std::vector<int> libMesh::ExodusII_IO_Helper::num_sides_per_set

Definition at line 455 of file exodusII_io_helper.h.

Referenced by read_sideset(), and read_sideset_info().

int libMesh::ExodusII_IO_Helper::num_time_steps

Definition at line 510 of file exodusII_io_helper.h.

Referenced by read_num_time_steps(), and read_time_steps().

bool libMesh::ExodusII_IO_Helper::opened_for_reading

Definition at line 550 of file exodusII_io_helper.h.

Referenced by close(), and open().

bool libMesh::ExodusII_IO_Helper::opened_for_writing
std::vector<int> libMesh::ExodusII_IO_Helper::side_list

Definition at line 470 of file exodusII_io_helper.h.

Referenced by read_sideset(), and read_sideset_info().

std::vector<int> libMesh::ExodusII_IO_Helper::ss_ids
std::vector<Real> libMesh::ExodusII_IO_Helper::time_steps

Definition at line 513 of file exodusII_io_helper.h.

Referenced by read_time_steps().

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

Definition at line 494 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 542 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_qa_records(), read_var_names_impl(), write_elements(), libMesh::Nemesis_IO_Helper::write_nodesets(), libMesh::Nemesis_IO_Helper::write_sidesets(), and write_var_names_impl().

std::vector<Real> libMesh::ExodusII_IO_Helper::x
std::vector<Real> libMesh::ExodusII_IO_Helper::y
std::vector<Real> libMesh::ExodusII_IO_Helper::z

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