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
John W. Peterson
Date
2002

Definition at line 86 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 575 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 1808 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().

1811 {
1812  // There may already be global variables in the file (for example,
1813  // if we're appending) and in that case, we
1814  // 1.) Cannot initialize them again.
1815  // 2.) Should check to be sure that the global variable names are the same.
1816 
1817  // Fills up names_from_file for us
1818  this->read_var_names(type);
1819 
1820  // Both the names of the global variables and their order must match
1821  if (names_from_file != names)
1822  {
1823  libMesh::err << "Error! The Exodus file already contains the variables:" << std::endl;
1824  for (std::size_t i=0; i<names_from_file.size(); ++i)
1825  libMesh::out << names_from_file[i] << std::endl;
1826 
1827  libMesh::err << "And you asked to write:" << std::endl;
1828  for (std::size_t i=0; i<names.size(); ++i)
1829  libMesh::out << names[i] << std::endl;
1830 
1831  libmesh_error_msg("Cannot overwrite existing variables in Exodus II file.");
1832  }
1833 }
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 787 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().

788 {
789  // Always call close on processor 0.
790  // If we're running on multiple processors, i.e. as one of several Nemesis files,
791  // we call close on all processors...
792  if ((this->processor_id() == 0) || (!_run_only_on_proc0))
793  {
794  // Don't close the file if it was never opened, this raises an Exodus error
796  {
798  EX_CHECK_ERR(ex_err, "Error closing Exodus file.");
799  message("Exodus file closed successfully.");
800  }
801  }
802 }
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::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::MeshTools::create_bounding_box(), libMesh::MeshTools::create_nodal_bounding_box(), libMesh::MeshRefinement::create_parent_error_vector(), libMesh::MeshTools::create_processor_bounding_box(), libMesh::MeshTools::create_subdomain_bounding_box(), 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::Nemesis_IO::read(), libMesh::CheckpointIO::read(), libMesh::XdrIO::read(), libMesh::CheckpointIO::read_header(), 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::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::MeshRefinement::uniformly_coarsen(), 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 1071 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.

1072 {
1073  // If we're processor 0, always create the file.
1074  // If we running on all procs, e.g. as one of several Nemesis files, also
1075  // call create there.
1076  if ((this->processor_id() == 0) || (!_run_only_on_proc0))
1077  {
1078  int
1079  comp_ws = 0,
1080  io_ws = 0;
1081 
1082  if (_single_precision)
1083  {
1084  comp_ws = cast_int<int>(sizeof(float));
1085  io_ws = cast_int<int>(sizeof(float));
1086  }
1087  // Fall back on double precision when necessary since ExodusII
1088  // doesn't seem to support long double
1089  else
1090  {
1091  comp_ws = cast_int<int>
1092  (std::min(sizeof(Real), sizeof(double)));
1093  io_ws = cast_int<int>
1094  (std::min(sizeof(Real), sizeof(double)));
1095  }
1096 
1097  ex_id = exII::ex_create(filename.c_str(), EX_CLOBBER, &comp_ws, &io_ws);
1098 
1099  EX_CHECK_ERR(ex_id, "Error creating ExodusII mesh file.");
1100 
1101  if (verbose)
1102  libMesh::out << "File created successfully." << std::endl;
1103  }
1104 
1105  opened_for_writing = true;
1106  current_filename = filename;
1107 }
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 536 of file exodusII_io_helper.C.

References block_ids.

Referenced by write_element_values().

537 {
538  libmesh_assert_less (static_cast<unsigned int>(index), block_ids.size());
539 
540  return block_ids[index];
541 }
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 545 of file exodusII_io_helper.C.

References block_ids, and id_to_block_names.

546 {
547  libmesh_assert_less (static_cast<unsigned int>(index), block_ids.size());
548 
549  return id_to_block_names[block_ids[index]];
550 }
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 2037 of file exodusII_io_helper.C.

2038 {
2039  std::vector<std::string>::const_iterator names_it = names.begin();
2040  std::vector<std::string>::const_iterator names_end = names.end();
2041 
2042  std::vector<std::string> complex_names;
2043 
2044  // This will loop over all names and create new "complex" names
2045  // (i.e. names that start with r_, i_ or a_
2046  for (; names_it != names_end; ++names_it)
2047  {
2048  std::stringstream name_real, name_imag, name_abs;
2049  name_real << "r_" << *names_it;
2050  name_imag << "i_" << *names_it;
2051  name_abs << "a_" << *names_it;
2052 
2053  complex_names.push_back(name_real.str());
2054  complex_names.push_back(name_imag.str());
2055  complex_names.push_back(name_abs.str());
2056  }
2057 
2058  return complex_names;
2059 }
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 572 of file exodusII_io_helper.C.

References nodeset_ids.

573 {
574  libmesh_assert_less (static_cast<unsigned int>(index), nodeset_ids.size());
575 
576  return nodeset_ids[index];
577 }
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 581 of file exodusII_io_helper.C.

References id_to_ns_names, and nodeset_ids.

582 {
583  libmesh_assert_less (static_cast<unsigned int>(index), nodeset_ids.size());
584 
585  return id_to_ns_names[nodeset_ids[index]];
586 }
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 554 of file exodusII_io_helper.C.

References ss_ids.

555 {
556  libmesh_assert_less (static_cast<unsigned int>(index), ss_ids.size());
557 
558  return ss_ids[index];
559 }
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 563 of file exodusII_io_helper.C.

References id_to_ss_names, and ss_ids.

564 {
565  libmesh_assert_less (static_cast<unsigned int>(index), ss_ids.size());
566 
567  return id_to_ss_names[ss_ids[index]];
568 }
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 1111 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::Elem::infinite(), 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().

1112 {
1113  // n_active_elem() is a parallel_only function
1114  unsigned int n_active_elem = mesh.n_active_elem();
1115 
1116  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1117  return;
1118 
1119  // If _write_as_dimension is nonzero, use it to set num_dim in the Exodus file.
1120  if (_write_as_dimension)
1123  num_dim = mesh.mesh_dimension();
1124  else
1125  num_dim = mesh.spatial_dimension();
1126 
1127  num_elem = mesh.n_elem();
1128 
1129  if (!use_discontinuous)
1130  {
1131  // Don't rely on mesh.n_nodes() here. If nodes have been
1132  // deleted recently, it will be incorrect.
1133  num_nodes = cast_int<int>(std::distance(mesh.nodes_begin(),
1134  mesh.nodes_end()));
1135  }
1136  else
1137  {
1140  for (; it!=end; ++it)
1141  num_nodes += (*it)->n_nodes();
1142  }
1143 
1144  std::vector<boundary_id_type> unique_side_boundaries;
1145  std::vector<boundary_id_type> unique_node_boundaries;
1146 
1147  mesh.get_boundary_info().build_side_boundary_ids(unique_side_boundaries);
1148  {
1149  // Add shell face boundaries to the list of side boundaries, since ExodusII
1150  // treats these the same way.
1151  std::vector<boundary_id_type> shellface_boundaries;
1152  mesh.get_boundary_info().build_shellface_boundary_ids(shellface_boundaries);
1153  for (std::size_t i=0; i<shellface_boundaries.size(); i++)
1154  unique_side_boundaries.push_back(shellface_boundaries[i]);
1155  }
1156  mesh.get_boundary_info().build_node_boundary_ids(unique_node_boundaries);
1157 
1158  num_side_sets = cast_int<int>(unique_side_boundaries.size());
1159  num_node_sets = cast_int<int>(unique_node_boundaries.size());
1160 
1161  //loop through element and map between block and element vector
1162  std::map<subdomain_id_type, std::vector<unsigned int> > subdomain_map;
1163 
1166  for (; it!=end; ++it)
1167  {
1168  const Elem * elem = *it;
1169 
1170  // We skip writing infinite elements to the Exodus file, so
1171  // don't put them in the subdomain_map. That way the number of
1172  // blocks should be correct.
1173 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1174  if (elem->infinite())
1175  continue;
1176 #endif
1177 
1178  subdomain_id_type cur_subdomain = elem->subdomain_id();
1179  subdomain_map[cur_subdomain].push_back(elem->id());
1180  }
1181  num_elem_blk = cast_int<int>(subdomain_map.size());
1182 
1183  if (str_title.size() > MAX_LINE_LENGTH)
1184  {
1185  libMesh::err << "Warning, Exodus files cannot have titles longer than "
1186  << MAX_LINE_LENGTH
1187  << " characters. Your title will be truncated."
1188  << std::endl;
1189  str_title.resize(MAX_LINE_LENGTH);
1190  }
1191 
1193  str_title.c_str(),
1194  num_dim,
1195  num_nodes,
1196  n_active_elem,
1197  num_elem_blk,
1198  num_node_sets,
1199  num_side_sets);
1200 
1201  EX_CHECK_ERR(ex_err, "Error initializing new Exodus file.");
1202 }
const BoundaryInfo & get_boundary_info() const
Definition: mesh_base.h:114
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 bool infinite() const =0
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:156
OStreamProxy err(std::cerr)
subdomain_id_type subdomain_id() const
Definition: elem.h:1805
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:147
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 1707 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().

1708 {
1709  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1710  return;
1711 
1712  // Quick return if there are no element 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_elem_vars > 0)
1722  {
1723  this->check_existing_vars(ELEMENTAL, names, this->elem_var_names);
1724  return;
1725  }
1726 
1727  // Set the flag so we can skip this stuff on subsequent calls to
1728  // initialize_element_variables()
1729  _elem_vars_initialized = true;
1730 
1731  this->write_var_names(ELEMENTAL, names);
1732 
1733  // Form the element variable truth table and send to Exodus.
1734  // This tells which variables are written to which blocks,
1735  // and can dramatically speed up writing element variables
1736  //
1737  // We really should initialize all entries in the truth table to 0
1738  // and then loop over all subdomains, setting their entries to 1
1739  // if a given variable exists on that subdomain. However,
1740  // we don't have that information, and the element variables
1741  // passed to us are padded with zeroes for the blocks where
1742  // they aren't defined. To be consistent with that, fill
1743  // the truth table with ones.
1744  std::vector<int> truth_tab(num_elem_blk*num_elem_vars, 1);
1746  num_elem_blk,
1747  num_elem_vars,
1748  &truth_tab[0]);
1749  EX_CHECK_ERR(ex_err, "Error writing element truth table.");
1750 }
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 1782 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().

1783 {
1784  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1785  return;
1786 
1787  // Quick return if there are no global variables to write
1788  if (names.size() == 0)
1789  return;
1790 
1792  return;
1793 
1794  // Be sure that variables in the file match what we are asking for
1795  if (num_global_vars > 0)
1796  {
1797  this->check_existing_vars(GLOBAL, names, this->global_var_names);
1798  return;
1799  }
1800 
1801  _global_vars_initialized = true;
1802 
1803  this->write_var_names(GLOBAL, names);
1804 }
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 1754 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().

1755 {
1756  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1757  return;
1758 
1759  // Quick return if there are no nodal variables to write
1760  if (names.size() == 0)
1761  return;
1762 
1763  // Quick return if we have already called this function
1765  return;
1766 
1767  // Be sure that variables in the file match what we are asking for
1768  if (num_nodal_vars > 0)
1769  {
1770  this->check_existing_vars(NODAL, names, this->nodal_var_names);
1771  return;
1772  }
1773 
1774  // Set the flag so we can skip the rest of this function on subsequent calls.
1775  _nodal_vars_initialized = true;
1776 
1777  this->write_var_names(NODAL, names);
1778 }
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 = "" 
)
Returns
The value obtained from a generic exII::ex_inquire() call.

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

807 {
808  int ret_int = 0;
809  char ret_char = 0;
810  float ret_float = 0.;
811 
813  req_info_in,
814  &ret_int,
815  &ret_float,
816  &ret_char);
817 
818  EX_CHECK_ERR(ex_err, error_msg);
819 
820  return ret_int;
821 }
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::MeshTools::create_processor_bounding_box(), 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::MeshRefinement::make_coarsening_compatible(), 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::NameBasedIO::read(), libMesh::Nemesis_IO::read(), libMesh::CheckpointIO::read(), libMesh::XdrIO::read(), libMesh::CheckpointIO::read_connectivity(), libMesh::CheckpointIO::read_nodes(), 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::MeshRefinement::uniformly_coarsen(), 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:722
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 _single_precision, 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  int comp_ws = 0;
336 
337  if (_single_precision)
338  comp_ws = cast_int<int>(sizeof(float));
339 
340  // Fall back on double precision when necessary since ExodusII
341  // doesn't seem to support long double
342  else
343  comp_ws = cast_int<int>(std::min(sizeof(Real), sizeof(double)));
344 
345  // Word size in bytes of the floating point data as they are stored
346  // in the ExodusII file. "If this argument is 0, the word size of the
347  // floating point data already stored in the file is returned"
348  int io_ws = 0;
349 
350  ex_id = exII::ex_open(filename,
351  read_only ? EX_READ : EX_WRITE,
352  &comp_ws,
353  &io_ws,
354  &ex_version);
355 
356  std::string err_msg = std::string("Error opening ExodusII mesh file: ") + std::string(filename);
357  EX_CHECK_ERR(ex_id, err_msg);
358  if (verbose) libMesh::out << "File opened successfully." << std::endl;
359 
360  if (read_only)
361  opened_for_reading = true;
362  else
363  opened_for_writing = true;
364 
365  current_filename = std::string(filename);
366 }
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 450 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.

451 {
452  if (verbose)
453  libMesh::out << "Title: \t" << &title[0] << std::endl
454  << "Mesh Dimension: \t" << num_dim << std::endl
455  << "Number of Nodes: \t" << num_nodes << std::endl
456  << "Number of elements: \t" << num_elem << std::endl
457  << "Number of elt blocks: \t" << num_elem_blk << std::endl
458  << "Number of node sets: \t" << num_node_sets << std::endl
459  << "Number of side sets: \t" << num_side_sets << std::endl;
460 }
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 501 of file exodusII_io_helper.C.

References num_nodes, x, y, and z.

502 {
503  for (int i=0; i<num_nodes; i++)
504  out_stream << "(" << x[i] << ", " << y[i] << ", " << z[i] << ")" << std::endl;
505 }
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::CheckpointIO::read(), libMesh::XdrIO::read(), read_elem_num_map(), libMesh::CheckpointIO::read_header(), 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::MeshRefinement::uniformly_coarsen(), libMesh::Parallel::Packing< Node * >::unpack(), libMesh::Parallel::Packing< Elem * >::unpack(), libMesh::DistributedMesh::update_parallel_id_counts(), libMesh::MeshTools::weight(), libMesh::NameBasedIO::write(), libMesh::CheckpointIO::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:720
void libMesh::ExodusII_IO_Helper::read_block_info ( )

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

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

510 {
511  block_ids.resize(num_elem_blk);
512  // Get all element block IDs.
514  block_ids.empty() ? libmesh_nullptr : &block_ids[0]);
515  // Usually, there is only one
516  // block since there is only
517  // one type of element.
518  // However, there could be more.
519 
520  EX_CHECK_ERR(ex_err, "Error getting block IDs.");
521  message("All block IDs retrieved successfully.");
522 
523  char name_buffer[MAX_STR_LENGTH+1];
524  for (int i=0; i<num_elem_blk; ++i)
525  {
527  block_ids[i], name_buffer);
528  EX_CHECK_ERR(ex_err, "Error getting block name.");
529  id_to_block_names[block_ids[i]] = name_buffer;
530  }
531  message("All block names retrieved successfully.");
532 }
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 591 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.

592 {
593  libmesh_assert_less (static_cast<unsigned int>(block), block_ids.size());
594 
596  block_ids[block],
597  &elem_type[0],
600  &num_attr);
601  if (verbose)
602  libMesh::out << "Reading a block of " << num_elem_this_blk
603  << " " << &elem_type[0] << "(s)"
604  << " having " << num_nodes_per_elem
605  << " nodes per element." << std::endl;
606 
607  EX_CHECK_ERR(ex_err, "Error getting block info.");
608  message("Info retrieved successfully for block: ", block);
609 
610 
611 
612  // Read in the connectivity of the elements of this block,
613  // watching out for the case where we actually have no
614  // elements in this block (possible with parallel files)
616 
617  if (!connect.empty())
618  {
619  ex_err = exII::ex_get_elem_conn(ex_id,
620  block_ids[block],
621  &connect[0]);
622 
623  EX_CHECK_ERR(ex_err, "Error reading block connectivity.");
624  message("Connectivity retrieved successfully for block: ", block);
625  }
626 }
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 631 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.

632 {
633  elem_num_map.resize(num_elem);
634 
636  elem_num_map.empty() ? libmesh_nullptr : &elem_num_map[0]);
637 
638  EX_CHECK_ERR(ex_err, "Error retrieving element number map.");
639  message("Element numbering map retrieved successfully.");
640 
641 
642  if (verbose)
643  {
644  libMesh::out << "[" << this->processor_id() << "] elem_num_map[i] = ";
645  for (unsigned int i=0; i<static_cast<unsigned int>(std::min(10, num_elem-1)); ++i)
646  libMesh::out << elem_num_map[i] << ", ";
647  libMesh::out << "... " << elem_num_map.back() << std::endl;
648  }
649 }
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 1004 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().

1007 {
1008  this->read_var_names(ELEMENTAL);
1009 
1010  // See if we can find the variable we are looking for
1011  std::size_t var_index = 0;
1012  bool found = false;
1013 
1014  // Do a linear search for nodal_var_name in nodal_var_names
1015  for (; var_index<elem_var_names.size(); ++var_index)
1016  if (elem_var_names[var_index] == elemental_var_name)
1017  {
1018  found = true;
1019  break;
1020  }
1021 
1022  if (!found)
1023  {
1024  libMesh::err << "Available variables: " << std::endl;
1025  for (std::size_t i=0; i<elem_var_names.size(); ++i)
1026  libMesh::err << elem_var_names[i] << std::endl;
1027 
1028  libmesh_error_msg("Unable to locate variable named: " << elemental_var_name);
1029  }
1030 
1031  // Sequential index which we can use to look up the element ID in the elem_num_map.
1032  unsigned ex_el_num = 0;
1033 
1034  for (unsigned i=0; i<static_cast<unsigned>(num_elem_blk); i++)
1035  {
1037  block_ids[i],
1041  libmesh_nullptr);
1042  EX_CHECK_ERR(ex_err, "Error getting number of elements in block.");
1043 
1044  std::vector<Real> block_elem_var_values(num_elem);
1046  time_step,
1047  var_index+1,
1048  block_ids[i],
1050  &block_elem_var_values[0]);
1051  EX_CHECK_ERR(ex_err, "Error getting elemental values.");
1052 
1053  for (unsigned j=0; j<static_cast<unsigned>(num_elem_this_blk); j++)
1054  {
1055  // Use the elem_num_map to obtain the ID of this element in the Exodus file,
1056  // and remember to subtract 1 since libmesh is zero-based and Exodus is 1-based.
1057  unsigned mapped_elem_id = this->elem_num_map[ex_el_num] - 1;
1058 
1059  // Store the elemental value in the map.
1060  elem_var_value_map[mapped_elem_id] = block_elem_var_values[j];
1061 
1062  // Go to the next sequential element ID.
1063  ex_el_num++;
1064  }
1065  }
1066 }
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 370 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.

371 {
373  &title[0],
374  &num_dim,
375  &num_nodes,
376  &num_elem,
377  &num_elem_blk,
378  &num_node_sets,
379  &num_side_sets);
380 
381  EX_CHECK_ERR(ex_err, "Error retrieving header info.");
382 
383  this->read_num_time_steps();
384 
386  EX_CHECK_ERR(ex_err, "Error reading number of nodal variables.");
387 
389  EX_CHECK_ERR(ex_err, "Error reading number of elemental variables.");
390 
392  EX_CHECK_ERR(ex_err, "Error reading number of global variables.");
393 
394  message("Exodus header info retrieved successfully.");
395 }
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 848 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().

849 {
850  // Read the nodal variable names from file, so we can see if we have the one we're looking for
851  this->read_var_names(NODAL);
852 
853  // See if we can find the variable we are looking for
854  std::size_t var_index = 0;
855  bool found = false;
856 
857  // Do a linear search for nodal_var_name in nodal_var_names
858  for (; var_index<nodal_var_names.size(); ++var_index)
859  {
860  found = (nodal_var_names[var_index] == nodal_var_name);
861  if (found)
862  break;
863  }
864 
865  if (!found)
866  {
867  libMesh::err << "Available variables: " << std::endl;
868  for (std::size_t i=0; i<nodal_var_names.size(); ++i)
869  libMesh::err << nodal_var_names[i] << std::endl;
870 
871  libmesh_error_msg("Unable to locate variable named: " << nodal_var_name);
872  }
873 
874  // Allocate enough space to store the nodal variable values
875  nodal_var_values.resize(num_nodes);
876 
877  // Call the Exodus API to read the nodal variable values
879  time_step,
880  var_index+1,
881  num_nodes,
882  &nodal_var_values[0]);
883  EX_CHECK_ERR(ex_err, "Error reading nodal variable values!");
884 }
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 481 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.

482 {
483  node_num_map.resize(num_nodes);
484 
486  node_num_map.empty() ? libmesh_nullptr : &node_num_map[0]);
487 
488  EX_CHECK_ERR(ex_err, "Error retrieving nodal number map.");
489  message("Nodal numbering map retrieved successfully.");
490 
491  if (verbose)
492  {
493  libMesh::out << "[" << this->processor_id() << "] node_num_map[i] = ";
494  for (unsigned int i=0; i<static_cast<unsigned int>(std::min(10, num_nodes-1)); ++i)
495  libMesh::out << node_num_map[i] << ", ";
496  libMesh::out << "... " << node_num_map.back() << std::endl;
497  }
498 }
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 464 of file exodusII_io_helper.C.

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

465 {
466  x.resize(num_nodes);
467  y.resize(num_nodes);
468  z.resize(num_nodes);
469 
471  static_cast<void *>(&x[0]),
472  static_cast<void *>(&y[0]),
473  static_cast<void *>(&z[0]));
474 
475  EX_CHECK_ERR(ex_err, "Error retrieving nodal data.");
476  message("Nodal data retrieved successfully.");
477 }
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 757 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.

758 {
759  libmesh_assert_less (static_cast<unsigned int>(id), nodeset_ids.size());
760  libmesh_assert_less (static_cast<unsigned int>(id), num_nodes_per_set.size());
761  libmesh_assert_less (static_cast<unsigned int>(id), num_node_df_per_set.size());
762 
764  nodeset_ids[id],
765  &num_nodes_per_set[id],
766  &num_node_df_per_set[id]);
767  EX_CHECK_ERR(ex_err, "Error retrieving nodeset parameters.");
768  message("Parameters retrieved successfully for nodeset: ", id);
769 
770  node_list.resize(num_nodes_per_set[id]);
771 
772  // Don't call ex_get_node_set unless there are actually nodes there to get.
773  // Exodus prints an annoying warning message in DEBUG mode otherwise...
774  if (num_nodes_per_set[id] > 0)
775  {
776  ex_err = exII::ex_get_node_set(ex_id,
777  nodeset_ids[id],
778  &node_list[0]);
779 
780  EX_CHECK_ERR(ex_err, "Error retrieving nodeset data.");
781  message("Data retrieved successfully for nodeset: ", id);
782  }
783 }
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 686 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.

687 {
688  nodeset_ids.resize(num_node_sets);
689  if (num_node_sets > 0)
690  {
692  &nodeset_ids[0]);
693  EX_CHECK_ERR(ex_err, "Error retrieving nodeset information.");
694  message("All nodeset information retrieved successfully.");
695 
696  // Resize appropriate data structures -- only do this once outnode the loop
699  }
700 
701  char name_buffer[MAX_STR_LENGTH+1];
702  for (int i=0; i<num_node_sets; ++i)
703  {
705  nodeset_ids[i], name_buffer);
706  EX_CHECK_ERR(ex_err, "Error getting node set name.");
707  id_to_ns_names[nodeset_ids[i]] = name_buffer;
708  }
709  message("All node set names retrieved successfully.");
710 }
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 840 of file exodusII_io_helper.C.

References EX_INQ_TIME, inquire(), and num_time_steps.

Referenced by read_header(), and read_time_steps().

841 {
843  this->inquire(exII::EX_INQ_TIME, "Error retrieving number of time steps");
844 }
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 400 of file exodusII_io_helper.C.

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

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

715 {
716  libmesh_assert_less (static_cast<unsigned int>(id), ss_ids.size());
717  libmesh_assert_less (static_cast<unsigned int>(id), num_sides_per_set.size());
718  libmesh_assert_less (static_cast<unsigned int>(id), num_df_per_set.size());
719  libmesh_assert_less_equal (static_cast<unsigned int>(offset), elem_list.size());
720  libmesh_assert_less_equal (static_cast<unsigned int>(offset), side_list.size());
721 
723  ss_ids[id],
724  &num_sides_per_set[id],
725  &num_df_per_set[id]);
726  EX_CHECK_ERR(ex_err, "Error retrieving sideset parameters.");
727  message("Parameters retrieved successfully for sideset: ", id);
728 
729 
730  // It's OK for offset==elem_list.size() as long as num_sides_per_set[id]==0
731  // because in that case we don't actually read anything...
732 #ifdef DEBUG
733  if (static_cast<unsigned int>(offset) == elem_list.size() ||
734  static_cast<unsigned int>(offset) == side_list.size() )
735  libmesh_assert_equal_to (num_sides_per_set[id], 0);
736 #endif
737 
738 
739  // Don't call ex_get_side_set unless there are actually sides there to get.
740  // Exodus prints an annoying warning in DEBUG mode otherwise...
741  if (num_sides_per_set[id] > 0)
742  {
743  ex_err = exII::ex_get_side_set(ex_id,
744  ss_ids[id],
745  &elem_list[offset],
746  &side_list[offset]);
747  EX_CHECK_ERR(ex_err, "Error retrieving sideset data.");
748  message("Data retrieved successfully for sideset: ", id);
749 
750  for (int i=0; i<num_sides_per_set[id]; i++)
751  id_list[i+offset] = ss_ids[id];
752  }
753 }
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 653 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.

654 {
655  ss_ids.resize(num_side_sets);
656  if (num_side_sets > 0)
657  {
659  &ss_ids[0]);
660  EX_CHECK_ERR(ex_err, "Error retrieving sideset information.");
661  message("All sideset information retrieved successfully.");
662 
663  // Resize appropriate data structures -- only do this once outside the loop
666 
667  // Inquire about the length of the concatenated side sets element list
668  num_elem_all_sidesets = inquire(exII::EX_INQ_SS_ELEM_LEN, "Error retrieving length of the concatenated side sets element list!");
669 
673  }
674 
675  char name_buffer[MAX_STR_LENGTH+1];
676  for (int i=0; i<num_side_sets; ++i)
677  {
679  ss_ids[i], name_buffer);
680  EX_CHECK_ERR(ex_err, "Error getting side set name.");
681  id_to_ss_names[ss_ids[i]] = name_buffer;
682  }
683  message("All side set names retrieved successfully.");
684 }
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 825 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.

826 {
827  // Make sure we have an up-to-date count of the number of time steps in the file.
828  this->read_num_time_steps();
829 
830  if (num_time_steps > 0)
831  {
832  time_steps.resize(num_time_steps);
834  EX_CHECK_ERR(ex_err, "Error reading timesteps!");
835  }
836 }
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 888 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().

889 {
890  switch (type)
891  {
892  case NODAL:
894  break;
895  case ELEMENTAL:
897  break;
898  case GLOBAL:
900  break;
901  default:
902  libmesh_error_msg("Unrecognized ExodusVarType " << type);
903  }
904 }
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 908 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().

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

References _coordinate_offset.

2032 {
2033  _coordinate_offset = p;
2034 }
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 2017 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 2024 of file exodusII_io_helper.C.

References _write_as_dimension.

2025 {
2026  _write_as_dimension = dim;
2027 }
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 1859 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::MeshBase::n_elem(), libMesh::MeshTools::n_elem(), num_elem_vars, libMesh::ParallelObject::processor_id(), and libMesh::Elem::subdomain_id().

1860 {
1861  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1862  return;
1863 
1864  // Loop over the element blocks and write the data one block at a time
1865  std::map<unsigned int, std::vector<unsigned int> > subdomain_map;
1866 
1867  // Ask the file how many element vars it has, store it in the num_elem_vars variable.
1869  EX_CHECK_ERR(ex_err, "Error reading number of elemental variables.");
1870 
1873 
1874  // loop through element and map between block and element vector
1875  for (; mesh_it!=end; ++mesh_it)
1876  {
1877  const Elem * elem = *mesh_it;
1878  subdomain_map[elem->subdomain_id()].push_back(elem->id());
1879  }
1880 
1881  // Use mesh.n_elem() to access into the values vector rather than
1882  // the number of elements the Exodus writer thinks the mesh has,
1883  // which may not include inactive elements.
1884  dof_id_type n_elem = mesh.n_elem();
1885 
1886  // For each variable, create a 'data' array which holds all the elemental variable
1887  // values *for a given block* on this processor, then write that data vector to file
1888  // before moving onto the next block.
1889  for (unsigned int i=0; i<static_cast<unsigned>(num_elem_vars); ++i)
1890  {
1891  // The size of the subdomain map is the number of blocks.
1892  std::map<unsigned int, std::vector<unsigned int> >::iterator it = subdomain_map.begin();
1893 
1894  for (unsigned int j=0; it!=subdomain_map.end(); ++it, ++j)
1895  {
1896  const std::vector<unsigned int> & elem_nums = (*it).second;
1897  const unsigned int num_elems_this_block =
1898  cast_int<unsigned int>(elem_nums.size());
1899  std::vector<Real> data(num_elems_this_block);
1900 
1901  for (unsigned int k=0; k<num_elems_this_block; ++k)
1902  data[k] = values[i*n_elem + elem_nums[k]];
1903 
1904  if (_single_precision)
1905  {
1906  std::vector<float> cast_data(data.begin(), data.end());
1907 
1909  timestep,
1910  i+1,
1911  this->get_block_id(j),
1912  num_elems_this_block,
1913  &cast_data[0]);
1914  }
1915  else
1916  {
1917  ex_err = exII::ex_put_elem_var(ex_id,
1918  timestep,
1919  i+1,
1920  this->get_block_id(j),
1921  num_elems_this_block,
1922  &data[0]);
1923  }
1924  EX_CHECK_ERR(ex_err, "Error writing element values.");
1925  }
1926  }
1927 
1929  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
1930 }
dof_id_type n_elem(const MeshBase::const_element_iterator &begin, const MeshBase::const_element_iterator &end)
Definition: mesh_tools.C:684
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:1805
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 1331 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_concat_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.

1332 {
1333  // n_active_elem() is a parallel_only function
1334  unsigned int n_active_elem = mesh.n_active_elem();
1335 
1336  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1337  return;
1338 
1339  // Map from block ID to a vector of element IDs in that block. Element
1340  // IDs are now of type dof_id_type, subdomain IDs are of type subdomain_id_type.
1341  typedef std::map<subdomain_id_type, std::vector<dof_id_type> > subdomain_map_type;
1342  subdomain_map_type subdomain_map;
1343 
1346  //loop through element and map between block and element vector
1347  for (; mesh_it!=end; ++mesh_it)
1348  {
1349  const Elem * elem = *mesh_it;
1350 
1351  // We skip writing infinite elements to the Exodus file, so
1352  // don't put them in the subdomain_map. That way the number of
1353  // blocks should be correct.
1354 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1355  if (elem->infinite())
1356  continue;
1357 #endif
1358 
1359  subdomain_map[ elem->subdomain_id() ].push_back(elem->id());
1360  }
1361 
1362  // element map vector
1363  num_elem_blk = cast_int<int>(subdomain_map.size());
1364  block_ids.resize(num_elem_blk);
1365  elem_num_map.resize(n_active_elem);
1366  std::vector<int>::iterator curr_elem_map_end = elem_num_map.begin();
1367 
1368  std::vector<int> elem_blk_id;
1369  std::vector<int> num_elem_this_blk_vec;
1370  std::vector<int> num_nodes_per_elem_vec;
1371  std::vector<int> num_attr_vec;
1372  NamesData elem_type_table(num_elem_blk, MAX_STR_LENGTH);
1373 
1374  // Note: It appears that there is a bug in exodusII::ex_put_name where
1375  // the index returned from the ex_id_lkup is erroneously used. For now
1376  // the work around is to use the alternative function ex_put_names, but
1377  // this function requires a char ** data structure.
1378  NamesData names_table(num_elem_blk, MAX_STR_LENGTH);
1379 
1380  // counter indexes into the block_ids vector
1381  unsigned int counter = 0;
1382  for (subdomain_map_type::iterator it=subdomain_map.begin(); it!=subdomain_map.end(); ++it, ++counter)
1383  {
1384  block_ids[counter] = (*it).first;
1385  names_table.push_back_entry(mesh.subdomain_name((*it).first));
1386 
1387  // Get a reference to a vector of element IDs for this subdomain.
1388  subdomain_map_type::mapped_type & tmp_vec = (*it).second;
1389 
1390  // Use the first element in this block to get representative information.
1391  // Note that Exodus assumes all elements in a block are of the same type!
1392  // We are using that same assumption here!
1394  const ExodusII_IO_Helper::Conversion conv =
1395  em.assign_conversion(mesh.elem_ref(tmp_vec[0]).type());
1396  num_nodes_per_elem = mesh.elem_ref(tmp_vec[0]).n_nodes();
1397 
1398  elem_blk_id.push_back((*it).first);
1399  elem_type_table.push_back_entry(conv.exodus_elem_type().c_str());
1400  num_elem_this_blk_vec.push_back(tmp_vec.size());
1401  num_nodes_per_elem_vec.push_back(num_nodes_per_elem);
1402  num_attr_vec.push_back(0); // we don't currently use elem block attributes.
1403  }
1404 
1405  // The "define_maps" parameter should be 0 if node_number_map and
1406  // elem_number_map will not be written later, and nonzero otherwise.
1408  &elem_blk_id[0],
1409  elem_type_table.get_char_star_star(),
1410  &num_elem_this_blk_vec[0],
1411  &num_nodes_per_elem_vec[0],
1412  &num_attr_vec[0],
1413  /*define_maps=*/0);
1414  EX_CHECK_ERR(ex_err, "Error writing element blocks.");
1415 
1416  // This counter is used to fill up the libmesh_elem_num_to_exodus map in the loop below.
1417  unsigned libmesh_elem_num_to_exodus_counter = 0;
1418 
1419  // node counter for discontinuous plotting
1420  unsigned int node_counter = 0;
1421  for (subdomain_map_type::iterator it=subdomain_map.begin(); it!=subdomain_map.end(); ++it)
1422  {
1423  // Get a reference to a vector of element IDs for this subdomain.
1424  subdomain_map_type::mapped_type & tmp_vec = (*it).second;
1425 
1426  //Use the first element in this block to get representative information.
1427  //Note that Exodus assumes all elements in a block are of the same type!
1428  //We are using that same assumption here!
1430  const ExodusII_IO_Helper::Conversion conv =
1431  em.assign_conversion(mesh.elem_ref(tmp_vec[0]).type());
1432  num_nodes_per_elem = mesh.elem_ref(tmp_vec[0]).n_nodes();
1433 
1434  connect.resize(tmp_vec.size()*num_nodes_per_elem);
1435 
1436  for (std::size_t i=0; i<tmp_vec.size(); i++)
1437  {
1438  unsigned int elem_id = tmp_vec[i];
1439  libmesh_elem_num_to_exodus[elem_id] = ++libmesh_elem_num_to_exodus_counter; // 1-based indexing for Exodus
1440 
1441  const Elem & elem = mesh.elem_ref(elem_id);
1442 
1443  // We *might* be able to get away with writing mixed element
1444  // types which happen to have the same number of nodes, but
1445  // do we actually *want* to get away with that?
1446  // .) No visualization software would be able to handle it.
1447  // .) There'd be no way for us to read it back in reliably.
1448  // .) Even elements with the same number of nodes may have different connectivities (?)
1449 
1450  // This needs to be more than an assert so we don't fail
1451  // with a mysterious segfault while trying to write mixed
1452  // element meshes in optimized mode.
1453  if (elem.type() != conv.get_canonical_type())
1454  libmesh_error_msg("Error: Exodus requires all elements with a given subdomain ID to be the same type.\n" \
1455  << "Can't write both " \
1456  << Utility::enum_to_string(elem.type()) \
1457  << " and " \
1459  << " in the same block!");
1460 
1461 
1462  for (unsigned int j=0; j<static_cast<unsigned int>(num_nodes_per_elem); ++j)
1463  {
1464  unsigned connect_index = (i*num_nodes_per_elem)+j;
1465  unsigned elem_node_index = conv.get_inverse_node_map(j); // inverse node map is for writing.
1466  if (verbose)
1467  {
1468  libMesh::out << "Exodus node index " << j
1469  << " = LibMesh node index " << elem_node_index << std::endl;
1470  }
1471 
1472  if (!use_discontinuous)
1473  {
1474  // The global id for the current node in libmesh.
1475  dof_id_type libmesh_node_id = elem.node_id(elem_node_index);
1476 
1477  // Find the zero-based libmesh id in the map, this
1478  // should be faster than doing linear searches on
1479  // the node_num_map.
1480  std::map<int, int>::iterator pos =
1481  libmesh_node_num_to_exodus.find(cast_int<int>(libmesh_node_id));
1482 
1483  // Make sure it was found.
1484  if (pos == libmesh_node_num_to_exodus.end())
1485  libmesh_error_msg("libmesh node id " << libmesh_node_id << " not found in node_num_map.");
1486 
1487  // Write the Exodus global node id associated with
1488  // this libmesh node number to the connectivity
1489  // array.
1490  connect[connect_index] = pos->second;
1491  }
1492  else
1493  {
1494  // FIXME: We are hard-coding the 1-based node
1495  // numbering assumption here, so writing
1496  // "discontinuous" Exodus files won't work with node
1497  // numberings that have "holes".
1498  connect[connect_index] = node_counter + elem_node_index + 1;
1499  }
1500  }
1501 
1502  node_counter += num_nodes_per_elem;
1503  }
1504 
1505  ex_err = exII::ex_put_elem_conn(ex_id, (*it).first, &connect[0]);
1506  EX_CHECK_ERR(ex_err, "Error writing element connectivities");
1507 
1508  // This transform command stores its result in a range that begins at the third argument,
1509  // so this command is adding values to the elem_num_map vector starting from curr_elem_map_end.
1510  curr_elem_map_end = std::transform(tmp_vec.begin(),
1511  tmp_vec.end(),
1512  curr_elem_map_end,
1513  std::bind2nd(std::plus<subdomain_map_type::mapped_type::value_type>(), 1)); // Adds one to each id to make a 1-based exodus file!
1514 
1515  // But if we don't want to add one, we just want to put the values
1516  // of tmp_vec into elem_map in the right location, we can use
1517  // std::copy().
1518  // curr_elem_map_end = std::copy(tmp_vec.begin(), tmp_vec.end(), curr_elem_map_end);
1519 
1520  counter++;
1521  }
1522 
1523  // write out the element number map that we created
1525  EX_CHECK_ERR(ex_err, "Error writing element map");
1526 
1527  // Write out the block names
1528  if (num_elem_blk > 0)
1529  {
1530  ex_err = exII::ex_put_names(ex_id, exII::EX_ELEM_BLOCK, names_table.get_char_star_star());
1531  EX_CHECK_ERR(ex_err, "Error writing element names");
1532  }
1533 
1534 }
virtual dof_id_type n_active_elem() const =0
virtual ElemType type() const =0
EXODUS_EXPORT int ex_put_concat_elem_block(int exoid, const void_int *elem_blk_id, char *elem_type[], const void_int *num_elem_this_blk, const void_int *num_nodes_per_elem, const void_int *num_attr, int define_maps)
The base class for all geometric element types.
Definition: elem.h:86
ExodusII_IO_Helper::Conversion assign_conversion(std::string type_str)
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:1805
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:483
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:1689
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 1995 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().

1996 {
1997  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1998  return;
1999 
2000  if (_single_precision)
2001  {
2002  std::vector<float> cast_values(values.begin(), values.end());
2003  ex_err = exII::ex_put_glob_vars(ex_id, timestep, num_global_vars, &cast_values[0]);
2004  }
2005  else
2006  {
2007  ex_err = exII::ex_put_glob_vars(ex_id, timestep, num_global_vars, &values[0]);
2008  }
2009  EX_CHECK_ERR(ex_err, "Error writing global values.");
2010 
2012  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
2013 }
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 1956 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().

1957 {
1958  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1959  return;
1960 
1961  // There may already be information records in the file (for
1962  // example, if we're appending) and in that case, according to the
1963  // Exodus documentation, writing more information records is not
1964  // supported.
1965  int num_info = inquire(exII::EX_INQ_INFO, "Error retrieving the number of information records from file!");
1966  if (num_info > 0)
1967  {
1968  libMesh::err << "Warning! The Exodus file already contains information records.\n"
1969  << "Exodus does not support writing additional records in this situation."
1970  << std::endl;
1971  return;
1972  }
1973 
1974  int num_records = cast_int<int>(records.size());
1975 
1976  if (num_records > 0)
1977  {
1978  NamesData info(num_records, MAX_LINE_LENGTH);
1979 
1980  // If an entry is longer than MAX_LINE_LENGTH characters it's not an error, we just
1981  // write the first MAX_LINE_LENGTH characters to the file.
1982  for (std::size_t i=0; i<records.size(); ++i)
1983  info.push_back_entry(records[i]);
1984 
1985  ex_err = exII::ex_put_info(ex_id, num_records, info.get_char_star_star());
1986  EX_CHECK_ERR(ex_err, "Error writing global values.");
1987 
1989  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
1990  }
1991 }
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 1206 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.

1207 {
1208  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1209  return;
1210 
1211  // Clear existing data from any previous calls.
1212  x.clear();
1213  y.clear();
1214  z.clear();
1215  node_num_map.clear();
1216 
1217  // Reserve space in the nodal coordinate vectors. num_nodes is
1218  // exact, this just allows us to do away with one potentially
1219  // error-inducing loop index.
1220  x.reserve(num_nodes);
1221  y.reserve(num_nodes);
1222  z.reserve(num_nodes);
1223 
1224  // And in the node_num_map - since the nodes aren't organized in
1225  // blocks, libmesh will always write out the identity map
1226  // here... unless there has been some refinement and coarsening, or
1227  // node deletion without a corresponding call to contract(). You
1228  // need to write this any time there could be 'holes' in the node
1229  // numbering, so we write it every time.
1230  node_num_map.reserve(num_nodes);
1231 
1232  // Clear out any previously-mapped node IDs.
1234 
1235  if (!use_discontinuous)
1236  {
1239  for (; it != end; ++it)
1240  {
1241  const Node & node = **it;
1242 
1243  x.push_back(node(0) + _coordinate_offset(0));
1244 
1245 #if LIBMESH_DIM > 1
1246  y.push_back(node(1) + _coordinate_offset(1));
1247 #else
1248  y.push_back(0.);
1249 #endif
1250 #if LIBMESH_DIM > 2
1251  z.push_back(node(2) + _coordinate_offset(2));
1252 #else
1253  z.push_back(0.);
1254 #endif
1255 
1256  // Fill in node_num_map entry with the proper (1-based) node id
1257  node_num_map.push_back(node.id() + 1);
1258 
1259  // Also map the zero-based libmesh node id to the 1-based
1260  // Exodus ID it will be assigned (this is equivalent to the
1261  // current size of the x vector).
1262  libmesh_node_num_to_exodus[ cast_int<int>(node.id()) ] = cast_int<int>(x.size());
1263  }
1264  }
1265  else
1266  {
1269 
1270  for (; it!=end; ++it)
1271  {
1272  const Elem * elem = *it;
1273 
1274  for (unsigned int n=0; n<elem->n_nodes(); n++)
1275  {
1276  x.push_back(elem->point(n)(0));
1277 #if LIBMESH_DIM > 1
1278  y.push_back(elem->point(n)(1));
1279 #else
1280  y.push_back(0.);
1281 #endif
1282 #if LIBMESH_DIM > 2
1283  z.push_back(elem->point(n)(2));
1284 #else
1285  z.push_back(0.);
1286 #endif
1287 
1288  // Let's skip the node_num_map in the discontinuous
1289  // case, since we're effectively duplicating nodes for
1290  // the sake of discontinuous visualization, so it isn't
1291  // clear how to deal with node_num_map here. This means
1292  // that writing discontinuous meshes won't work with
1293  // element numberings that have "holes".
1294  }
1295  }
1296  }
1297 
1298  if (_single_precision)
1299  {
1300  std::vector<float>
1301  x_single(x.begin(), x.end()),
1302  y_single(y.begin(), y.end()),
1303  z_single(z.begin(), z.end());
1304 
1306  x_single.empty() ? libmesh_nullptr : &x_single[0],
1307  y_single.empty() ? libmesh_nullptr : &y_single[0],
1308  z_single.empty() ? libmesh_nullptr : &z_single[0]);
1309  }
1310  else
1311  {
1312  ex_err = exII::ex_put_coord(ex_id,
1313  x.empty() ? libmesh_nullptr : &x[0],
1314  y.empty() ? libmesh_nullptr : &y[0],
1315  z.empty() ? libmesh_nullptr : &z[0]);
1316  }
1317 
1318 
1319  EX_CHECK_ERR(ex_err, "Error writing coordinates to Exodus file.");
1320 
1321  if (!use_discontinuous)
1322  {
1323  // Also write the (1-based) node_num_map to the file.
1325  EX_CHECK_ERR(ex_err, "Error writing node_num_map");
1326  }
1327 }
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:1667
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 1934 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().

1935 {
1936  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1937  return;
1938 
1939  if (_single_precision)
1940  {
1941  std::vector<float> cast_values(values.begin(), values.end());
1942  ex_err = exII::ex_put_nodal_var(ex_id, timestep, var_id, num_nodes, &cast_values[0]);
1943  }
1944  else
1945  {
1946  ex_err = exII::ex_put_nodal_var(ex_id, timestep, var_id, num_nodes, &values[0]);
1947  }
1948  EX_CHECK_ERR(ex_err, "Error writing nodal values.");
1949 
1951  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
1952 }
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 1659 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().

1660 {
1661  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1662  return;
1663 
1664  std::vector< dof_id_type > nl;
1665  std::vector< boundary_id_type > il;
1666 
1667  mesh.get_boundary_info().build_node_list(nl, il);
1668 
1669  // Maps from nodeset id to the nodes
1670  std::map<boundary_id_type, std::vector<int> > node;
1671 
1672  // Accumulate the vectors to pass into ex_put_node_set
1673  for (std::size_t i=0; i<nl.size(); i++)
1674  node[il[i]].push_back(nl[i]+1);
1675 
1676  std::vector<boundary_id_type> node_boundary_ids;
1677  mesh.get_boundary_info().build_node_boundary_ids(node_boundary_ids);
1678 
1679  // Write out the nodeset names, but only if there is something to write
1680  if (node_boundary_ids.size() > 0)
1681  {
1682  NamesData names_table(node_boundary_ids.size(), MAX_STR_LENGTH);
1683 
1684  for (std::size_t i=0; i<node_boundary_ids.size(); i++)
1685  {
1686  boundary_id_type nodeset_id = node_boundary_ids[i];
1687 
1688  int actual_id = nodeset_id;
1689 
1690  names_table.push_back_entry(mesh.get_boundary_info().get_nodeset_name(nodeset_id));
1691 
1692  ex_err = exII::ex_put_node_set_param(ex_id, actual_id, node[nodeset_id].size(), 0);
1693  EX_CHECK_ERR(ex_err, "Error writing nodeset parameters");
1694 
1695  ex_err = exII::ex_put_node_set(ex_id, actual_id, &node[nodeset_id][0]);
1696  EX_CHECK_ERR(ex_err, "Error writing nodesets");
1697  }
1698 
1699  // Write out the nodeset names
1700  ex_err = exII::ex_put_names(ex_id, exII::EX_NODE_SET, names_table.get_char_star_star());
1701  EX_CHECK_ERR(ex_err, "Error writing nodeset names");
1702  }
1703 }
const BoundaryInfo & get_boundary_info() const
Definition: mesh_base.h:114
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 1539 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_sets(), 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_nullptr, libMesh::ParallelObject::processor_id(), libMesh::ExodusII_IO_Helper::NamesData::push_back_entry(), and side.

1540 {
1541  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1542  return;
1543 
1545 
1546  // Maps from sideset id to the element and sides
1547  std::map<int, std::vector<int> > elem;
1548  std::map<int, std::vector<int> > side;
1549  std::vector<boundary_id_type> side_boundary_ids;
1550 
1551  {
1552  std::vector< dof_id_type > el;
1553  std::vector< unsigned short int > sl;
1554  std::vector< boundary_id_type > il;
1555 
1556  mesh.get_boundary_info().build_side_list(el, sl, il);
1557 
1558  // Accumulate the vectors to pass into ex_put_side_set
1559  for (std::size_t i=0; i<el.size(); i++)
1560  {
1561  std::vector<const Elem *> family;
1562 #ifdef LIBMESH_ENABLE_AMR
1563 
1567  mesh.elem_ref(el[i]).active_family_tree_by_side(family, sl[i], false);
1568 #else
1569  family.push_back(mesh.elem_ptr(el[i]));
1570 #endif
1571 
1572  for (std::size_t j=0; j<family.size(); ++j)
1573  {
1574  const ExodusII_IO_Helper::Conversion conv =
1575  em.assign_conversion(mesh.elem_ptr(family[j]->id())->type());
1576 
1577  // Use the libmesh to exodus data structure map to get the proper sideset IDs
1578  // The data structure contains the "collapsed" contiguous ids
1579  elem[il[i]].push_back(libmesh_elem_num_to_exodus[family[j]->id()]);
1580  side[il[i]].push_back(conv.get_inverse_side_map(sl[i]));
1581  }
1582  }
1583 
1584  mesh.get_boundary_info().build_side_boundary_ids(side_boundary_ids);
1585  }
1586 
1587  {
1588  // add data for shell faces, if needed
1589 
1590  std::vector< dof_id_type > el;
1591  std::vector< unsigned short int > sl;
1592  std::vector< boundary_id_type > il;
1593 
1594  mesh.get_boundary_info().build_shellface_list(el, sl, il);
1595 
1596  // Accumulate the vectors to pass into ex_put_side_set
1597  for (std::size_t i=0; i<el.size(); i++)
1598  {
1599  std::vector<const Elem *> family;
1600 #ifdef LIBMESH_ENABLE_AMR
1601 
1605  mesh.elem_ref(el[i]).active_family_tree_by_side(family, sl[i], false);
1606 #else
1607  family.push_back(mesh.elem_ptr(el[i]));
1608 #endif
1609 
1610  for (std::size_t j=0; j<family.size(); ++j)
1611  {
1612  const ExodusII_IO_Helper::Conversion conv =
1613  em.assign_conversion(mesh.elem_ptr(family[j]->id())->type());
1614 
1615  // Use the libmesh to exodus data structure map to get the proper sideset IDs
1616  // The data structure contains the "collapsed" contiguous ids
1617  elem[il[i]].push_back(libmesh_elem_num_to_exodus[family[j]->id()]);
1618  side[il[i]].push_back(conv.get_inverse_shellface_map(sl[i]));
1619  }
1620  }
1621 
1622  std::vector<boundary_id_type> shellface_boundary_ids;
1623  mesh.get_boundary_info().build_shellface_boundary_ids(shellface_boundary_ids);
1624  for (std::size_t i=0; i<shellface_boundary_ids.size(); i++)
1625  side_boundary_ids.push_back(shellface_boundary_ids[i]);
1626  }
1627 
1628  // Write out the sideset names, but only if there is something to write
1629  if (side_boundary_ids.size() > 0)
1630  {
1631  NamesData names_table(side_boundary_ids.size(), MAX_STR_LENGTH);
1632 
1633  std::vector<exII::ex_set> sets(side_boundary_ids.size());
1634 
1635  for (std::size_t i=0; i<side_boundary_ids.size(); i++)
1636  {
1637  boundary_id_type ss_id = side_boundary_ids[i];
1638  names_table.push_back_entry(mesh.get_boundary_info().get_sideset_name(ss_id));
1639 
1640  sets[i].id = ss_id;
1641  sets[i].type = exII::EX_SIDE_SET;
1642  sets[i].num_entry = elem[ss_id].size();
1643  sets[i].num_distribution_factor = 0;
1644  sets[i].entry_list = &elem[ss_id][0];
1645  sets[i].extra_list = &side[ss_id][0];
1646  sets[i].distribution_factor_list = libmesh_nullptr;
1647  }
1648 
1649  ex_err = exII::ex_put_sets(ex_id, side_boundary_ids.size(), &sets[0]);
1650  EX_CHECK_ERR(ex_err, "Error writing sidesets");
1651 
1652  ex_err = exII::ex_put_names(ex_id, exII::EX_SIDE_SET, names_table.get_char_star_star());
1653  EX_CHECK_ERR(ex_err, "Error writing sideset names");
1654  }
1655 }
const BoundaryInfo & get_boundary_info() const
Definition: mesh_base.h:114
unsigned short int side
Definition: xdr_io.C:49
const class libmesh_nullptr_t libmesh_nullptr
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
EXODUS_EXPORT int ex_put_sets(int exoid, size_t set_count, const struct ex_set *sets)
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_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:483
void active_family_tree_by_side(std::vector< const Elem * > &family, const unsigned int side, const bool reset=true) const
Definition: elem.C:1767
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 1837 of file exodusII_io_helper.C.

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

1838 {
1839  if ((_run_only_on_proc0) && (this->processor_id() != 0))
1840  return;
1841 
1842  if (_single_precision)
1843  {
1844  float cast_time = time;
1845  ex_err = exII::ex_put_time(ex_id, timestep, &cast_time);
1846  }
1847  else
1848  {
1849  ex_err = exII::ex_put_time(ex_id, timestep, &time);
1850  }
1851  EX_CHECK_ERR(ex_err, "Error writing timestep.");
1852 
1854  EX_CHECK_ERR(ex_err, "Error flushing buffers to file.");
1855 }
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 948 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().

949 {
950  switch (type)
951  {
952  case NODAL:
953  this->write_var_names_impl("n", num_nodal_vars, names);
954  break;
955  case ELEMENTAL:
956  this->write_var_names_impl("e", num_elem_vars, names);
957  break;
958  case GLOBAL:
959  this->write_var_names_impl("g", num_global_vars, names);
960  break;
961  default:
962  libmesh_error_msg("Unrecognized ExodusVarType " << type);
963  }
964 }
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 968 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().

969 {
970  // Update the count variable so that it's available to other parts of the class.
971  count = cast_int<int>(names.size());
972 
973  // Write that number of variables to the file.
974  ex_err = exII::ex_put_var_param(ex_id, var_type, count);
975  EX_CHECK_ERR(ex_err, "Error setting number of vars.");
976 
977  if (names.size() > 0)
978  {
979  NamesData names_table(names.size(), MAX_STR_LENGTH);
980 
981  // Store the input names in the format required by Exodus.
982  for (std::size_t i=0; i<names.size(); ++i)
983  names_table.push_back_entry(names[i]);
984 
985  if (verbose)
986  {
987  libMesh::out << "Writing variable name(s) to file: " << std::endl;
988  for (std::size_t i=0; i<names.size(); ++i)
989  libMesh::out << names_table.get_char_star(i) << std::endl;
990  }
991 
992  ex_err = exII::ex_put_var_names(ex_id,
993  var_type,
994  cast_int<int>(names.size()),
995  names_table.get_char_star_star()
996  );
997 
998  EX_CHECK_ERR(ex_err, "Error writing variable names.");
999  }
1000 }
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 601 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 583 of file exodusII_io_helper.h.

Referenced by initialize_element_variables().

bool libMesh::ExodusII_IO_Helper::_global_vars_initialized
protected

Definition at line 586 of file exodusII_io_helper.h.

Referenced by initialize_global_variables().

bool libMesh::ExodusII_IO_Helper::_nodal_vars_initialized
protected

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

Referenced by initialize(), and write_as_dimension().

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

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

Referenced by create(), and open().

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

Definition at line 472 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 502 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 536 of file exodusII_io_helper.h.

int libMesh::ExodusII_IO_Helper::ex_id

Definition at line 409 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 539 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 481 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 542 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 544 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 543 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 527 of file exodusII_io_helper.h.

Referenced by read_nodal_var_values().

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

Definition at line 478 of file exodusII_io_helper.h.

Referenced by read_nodeset().

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

Definition at line 484 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 442 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 466 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 445 of file exodusII_io_helper.h.

Referenced by read_sideset_info().

int libMesh::ExodusII_IO_Helper::num_elem_this_blk

Definition at line 436 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 469 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 463 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 460 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 515 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 555 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 475 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 518 of file exodusII_io_helper.h.

Referenced by read_time_steps().

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

Definition at line 499 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 547 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: