libMesh::Nemesis_IO Class Reference
#include <nemesis_io.h>

Public Member Functions | |
| Nemesis_IO (MeshBase &mesh) | |
| virtual | ~Nemesis_IO () |
| virtual void | read (const std::string &base_filename) |
| virtual void | write (const std::string &base_filename) |
| void | write_timestep (const std::string &fname, const EquationSystems &es, const int timestep, const Real time) |
| void | write_nodal_data (const std::string &fname, const std::vector< Number > &soln, const std::vector< std::string > &names) |
| void | verbose (bool set_verbosity) |
| void | write_global_data (const std::vector< Number > &, const std::vector< std::string > &) |
| void | write_information_records (const std::vector< std::string > &) |
| virtual void | write_equation_systems (const std::string &, const EquationSystems &, const std::set< std::string > *system_names=NULL) |
| unsigned int & | ascii_precision () |
Protected Member Functions | |
| MeshBase & | mesh () |
| void | skip_comment_lines (std::istream &in, const char comment_start) |
| const MeshBase & | mesh () const |
Protected Attributes | |
| std::vector< bool > | elems_of_dimension |
| const bool | _is_parallel_format |
Private Attributes | |
| Nemesis_IO_Helper * | nemhelper |
| int | _timestep |
| bool | _verbose |
Detailed Description
The Nemesis_IO class implements reading parallel meshes in the Nemesis file format from Sandia National Labs. Nemesis files are essentially in the Exodus format plus some additional information. All the Nemesis files for a single mesh have the same basename, e.g. cylinder.e, followed by ".size.rank", where size is the total number of files the Mesh is split into and rank is the ID of the processor's elements that were written to the file.
Definition at line 52 of file nemesis_io.h.
Constructor & Destructor Documentation
| libMesh::Nemesis_IO::Nemesis_IO | ( | MeshBase & | mesh | ) | [explicit] |
Constructor. Takes a writeable reference to a mesh object. This is the constructor required to read a mesh.
Definition at line 80 of file nemesis_io.C.
00080 : 00081 MeshInput<MeshBase> (mesh, /*is_parallel_format=*/true), 00082 MeshOutput<MeshBase> (mesh, /*is_parallel_format=*/true), 00083 #if defined(LIBMESH_HAVE_EXODUS_API) && defined(LIBMESH_HAVE_NEMESIS_API) 00084 nemhelper(new Nemesis_IO_Helper), 00085 #endif 00086 _timestep(1), 00087 _verbose (false) 00088 { 00089 }
| libMesh::Nemesis_IO::~Nemesis_IO | ( | ) | [virtual] |
Destructor.
Definition at line 92 of file nemesis_io.C.
References nemhelper.
00093 { 00094 #if defined(LIBMESH_HAVE_EXODUS_API) && defined(LIBMESH_HAVE_NEMESIS_API) 00095 00096 delete nemhelper; 00097 #endif 00098 }
Member Function Documentation
| unsigned int& libMesh::MeshOutput< MeshBase >::ascii_precision | ( | ) | [inherited] |
Return/set the precision to use when writing ASCII files.
By default we use numeric_limits<Real>::digits10 + 2, which should be enough to write out to ASCII and get the exact same Real back when reading in.
Referenced by libMesh::TecplotIO::write_ascii(), libMesh::GMVIO::write_ascii_new_impl(), and libMesh::GMVIO::write_ascii_old_impl().
| const MeshBase & libMesh::MeshOutput< MeshBase >::mesh | ( | ) | const [protected, inherited] |
Returns the object as a read-only reference.
Referenced by libMesh::PostscriptIO::write(), libMesh::FroIO::write(), libMesh::EnsightIO::write(), libMesh::DivaIO::write(), libMesh::TecplotIO::write_ascii(), libMesh::MEDITIO::write_ascii(), libMesh::TecplotIO::write_binary(), libMesh::EnsightIO::write_geometry_ascii(), libMesh::EnsightIO::write_scalar_ascii(), libMesh::GnuPlotIO::write_solution(), libMesh::DivaIO::write_stream(), and libMesh::EnsightIO::write_vector_ascii().
| MeshBase & libMesh::MeshInput< MeshBase >::mesh | ( | ) | [protected, inherited] |
Returns the object as a writeable reference.
Referenced by libMesh::GMVIO::_read_materials(), libMesh::GMVIO::_read_nodes(), libMesh::GMVIO::_read_one_cell(), libMesh::AbaqusIO::assign_boundary_node_ids(), libMesh::AbaqusIO::assign_sideset_ids(), libMesh::AbaqusIO::assign_subdomain_ids(), libMesh::VTKIO::cells_to_vtk(), libMesh::GMVIO::copy_nodal_solution(), libMesh::UNVIO::element_in(), libMesh::TetGenIO::element_in(), libMesh::UNVIO::element_out(), libMesh::UNVIO::node_in(), libMesh::TetGenIO::node_in(), libMesh::UNVIO::node_out(), libMesh::VTKIO::nodes_to_vtk(), libMesh::XdrIO::read(), libMesh::VTKIO::read(), libMesh::TetGenIO::read(), read(), libMesh::GMVIO::read(), libMesh::ExodusII_IO::read(), libMesh::AbaqusIO::read(), libMesh::LegacyXdrIO::read_ascii(), libMesh::AbaqusIO::read_elements(), libMesh::UNVIO::read_implementation(), libMesh::UCDIO::read_implementation(), libMesh::LegacyXdrIO::read_mesh(), libMesh::GmshIO::read_mesh(), libMesh::AbaqusIO::read_nodes(), libMesh::XdrIO::read_serialized_bcs(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::XdrIO::read_serialized_nodes(), libMesh::OFFIO::read_stream(), libMesh::MatlabIO::read_stream(), libMesh::XdrIO::write(), libMesh::TetGenIO::write(), write(), libMesh::GMVIO::write_ascii_new_impl(), libMesh::GMVIO::write_ascii_old_impl(), libMesh::GMVIO::write_binary(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::UNVIO::write_implementation(), libMesh::UCDIO::write_implementation(), libMesh::LegacyXdrIO::write_mesh(), libMesh::GmshIO::write_mesh(), libMesh::UCDIO::write_nodal_data(), write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data_discontinuous(), libMesh::XdrIO::write_parallel(), libMesh::GmshIO::write_post(), libMesh::XdrIO::write_serialized_bcs(), libMesh::XdrIO::write_serialized_nodes(), and libMesh::LegacyXdrIO::write_soln().
| void libMesh::Nemesis_IO::read | ( | const std::string & | base_filename | ) | [virtual] |
Implements reading the mesh from several different files. You provide the basename, then LibMesh appends the ".size.rank" depending on libMesh::n_processors() and libMesh::processor_id().
Implements libMesh::MeshInput< MeshBase >.
Definition at line 116 of file nemesis_io.C.
References _verbose, libMesh::MeshBase::add_elem(), libMesh::MeshBase::add_point(), libMesh::Parallel::Communicator::allgather(), libMesh::Parallel::Communicator::alltoall(), libMesh::Parallel::any_source, libMesh::ExodusII_IO_Helper::ElementMaps::assign_conversion(), libMesh::ExodusII_IO_Helper::block_ids, libMesh::MeshBase::boundary_info, libMesh::Elem::build(), libMesh::Parallel::Communicator_World, libMesh::CommWorld, libMesh::ExodusII_IO_Helper::connect, libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::Utility::deallocate(), libMesh::MeshBase::delete_remote_elements(), libMesh::Elem::dim(), libMesh::MeshBase::elem(), libMesh::ExodusII_IO_Helper::elem_type, libMesh::MeshInput< MeshBase >::elems_of_dimension, libMesh::err, libMesh::ExodusII_IO_Helper::Conversion::get_canonical_type(), libMesh::Nemesis_IO_Helper::get_cmap_params(), libMesh::Nemesis_IO_Helper::get_eb_info_global(), libMesh::Nemesis_IO_Helper::get_elem_cmap(), libMesh::ExodusII_IO_Helper::get_elem_list(), libMesh::ExodusII_IO_Helper::get_id_list(), libMesh::Nemesis_IO_Helper::get_init_global(), libMesh::Nemesis_IO_Helper::get_loadbal_param(), libMesh::Nemesis_IO_Helper::get_node_cmap(), libMesh::ExodusII_IO_Helper::get_node_list(), libMesh::ExodusII_IO_Helper::Conversion::get_node_map(), libMesh::Nemesis_IO_Helper::get_node_map(), libMesh::ExodusII_IO_Helper::get_nodeset_id(), libMesh::Nemesis_IO_Helper::get_ns_param_global(), libMesh::ExodusII_IO_Helper::get_num_node_sets(), libMesh::ExodusII_IO_Helper::get_num_side_sets(), libMesh::ExodusII_IO_Helper::get_num_sides_per_set(), libMesh::ExodusII_IO_Helper::get_side_list(), libMesh::ExodusII_IO_Helper::Conversion::get_side_map(), libMesh::Nemesis_IO_Helper::get_ss_param_global(), libMesh::Nemesis_IO_Helper::global_sideset_ids, libMesh::DofObject::id(), libMesh::MeshTools::Generation::Private::idx(), libMesh::MeshBase::is_serial(), libMesh::Parallel::Communicator::max(), libMesh::MeshInput< MeshBase >::mesh(), libMesh::MeshBase::mesh_dimension(), std::min(), libMesh::MeshBase::n_elem(), libMesh::MeshBase::n_nodes(), libMesh::n_processors(), nemhelper, libMesh::Nemesis_IO_Helper::node_cmap_ids, libMesh::Nemesis_IO_Helper::node_cmap_node_cnts, libMesh::Nemesis_IO_Helper::node_cmap_node_ids, libMesh::Nemesis_IO_Helper::node_cmap_proc_ids, libMesh::Nemesis_IO_Helper::node_mapb, libMesh::Nemesis_IO_Helper::node_mape, libMesh::Nemesis_IO_Helper::node_mapi, libMesh::ExodusII_IO_Helper::node_num_map, libMesh::MeshBase::node_ptr(), libMesh::Nemesis_IO_Helper::num_border_elems, libMesh::Nemesis_IO_Helper::num_border_nodes, libMesh::ExodusII_IO_Helper::num_elem, libMesh::ExodusII_IO_Helper::num_elem_all_sidesets, libMesh::ExodusII_IO_Helper::num_elem_blk, libMesh::Nemesis_IO_Helper::num_elem_cmaps, libMesh::ExodusII_IO_Helper::num_elem_this_blk, libMesh::Nemesis_IO_Helper::num_elems_global, libMesh::Nemesis_IO_Helper::num_external_nodes, libMesh::Nemesis_IO_Helper::num_global_side_counts, libMesh::Nemesis_IO_Helper::num_internal_elems, libMesh::Nemesis_IO_Helper::num_internal_nodes, libMesh::Nemesis_IO_Helper::num_node_cmaps, libMesh::ExodusII_IO_Helper::num_node_sets, libMesh::ExodusII_IO_Helper::num_nodes, libMesh::Nemesis_IO_Helper::num_nodes_global, libMesh::ExodusII_IO_Helper::num_nodes_per_elem, libMesh::ExodusII_IO_Helper::num_side_sets, libMesh::ExodusII_IO_Helper::open(), libMesh::out, libMesh::MeshBase::parallel_n_elem(), libMesh::MeshBase::parallel_n_nodes(), libMesh::ExodusII_IO_Helper::print_header(), libMesh::Parallel::Communicator::probe(), libMesh::DofObject::processor_id(), libMesh::processor_id(), libMesh::ExodusII_IO_Helper::read_block_info(), libMesh::ExodusII_IO_Helper::read_elem_in_block(), libMesh::ExodusII_IO_Helper::read_elem_num_map(), libMesh::ExodusII_IO_Helper::read_header(), libMesh::ExodusII_IO_Helper::read_node_num_map(), libMesh::ExodusII_IO_Helper::read_nodes(), libMesh::ExodusII_IO_Helper::read_nodeset(), libMesh::ExodusII_IO_Helper::read_nodeset_info(), libMesh::ExodusII_IO_Helper::read_sideset(), libMesh::ExodusII_IO_Helper::read_sideset_info(), libMesh::Parallel::Communicator::receive(), libMesh::Parallel::Communicator::send(), libMesh::DofObject::set_id(), libMesh::MeshBase::set_mesh_dimension(), libMesh::Elem::set_node(), libMesh::Elem::subdomain_id(), libMesh::Parallel::Communicator::sum(), libMesh::Elem::type(), libMesh::MeshBase::update_post_partitioning(), libMesh::Parallel::wait(), libMesh::ExodusII_IO_Helper::x, libMesh::ExodusII_IO_Helper::y, and libMesh::ExodusII_IO_Helper::z.
00117 { 00118 // On one processor, Nemesis and ExodusII should be equivalent, so 00119 // let's cowardly defer to that implementation... 00120 if (libMesh::n_processors() == 1) 00121 { 00122 // We can do this in one line but if the verbose flag was set in this 00123 // object, it will no longer be set... thus no extra print-outs for serial runs. 00124 // ExodusII_IO(this->mesh()).read (base_filename); // ambiguous when Nemesis_IO is multiply-inherited 00125 00126 MeshBase& mesh = MeshInput<MeshBase>::mesh(); 00127 ExodusII_IO(mesh).read (base_filename); 00128 return; 00129 } 00130 00131 START_LOG ("read()","Nemesis_IO"); 00132 00133 // This function must be run on all processors at once 00134 parallel_only(); 00135 00136 if (_verbose) 00137 { 00138 libMesh::out << "[" << libMesh::processor_id() << "] "; 00139 libMesh::out << "Reading Nemesis file on processor: " << libMesh::processor_id() << std::endl; 00140 } 00141 00142 // Construct the Nemesis filename based on the number of processors and the 00143 // current processor ID. 00144 std::string nemesis_filename = nemhelper->construct_nemesis_filename(base_filename); 00145 00146 if (_verbose) 00147 libMesh::out << "Opening file: " << nemesis_filename << std::endl; 00148 00149 // Open the Exodus file 00150 nemhelper->open(nemesis_filename.c_str()); 00151 00152 // Get a reference to the mesh. We need to be specific 00153 // since Nemesis_IO is multiply-inherited 00154 // MeshBase& mesh = this->mesh(); 00155 MeshBase& mesh = MeshInput<MeshBase>::mesh(); 00156 00157 // Local information: Read the following information from the standard Exodus header 00158 // title[0] 00159 // num_dim 00160 // num_nodes 00161 // num_elem 00162 // num_elem_blk 00163 // num_node_sets 00164 // num_side_sets 00165 nemhelper->read_header(); 00166 nemhelper->print_header(); 00167 00168 // Get global information: number of nodes, elems, blocks, nodesets and sidesets 00169 nemhelper->get_init_global(); 00170 00171 // Get "load balance" information. This includes the number of internal & border 00172 // nodes and elements as well as the number of communication maps. 00173 nemhelper->get_loadbal_param(); 00174 00175 // Do some error checking 00176 if (nemhelper->num_external_nodes) 00177 { 00178 libMesh::err << "ERROR: there should be no external nodes in an element-based partitioning!" 00179 << std::endl; 00180 libmesh_error(); 00181 } 00182 00183 libmesh_assert_equal_to (nemhelper->num_nodes, 00184 (nemhelper->num_internal_nodes + 00185 nemhelper->num_border_nodes)); 00186 00187 libmesh_assert_equal_to (nemhelper->num_elem, 00188 (nemhelper->num_internal_elems + 00189 nemhelper->num_border_elems)); 00190 00191 libmesh_assert_less_equal (nemhelper->num_nodes, nemhelper->num_nodes_global); 00192 libmesh_assert_less_equal (nemhelper->num_elem, nemhelper->num_elems_global); 00193 00194 // Read nodes from the exodus file: this fills the nemhelper->x,y,z arrays. 00195 nemhelper->read_nodes(); 00196 00197 // Reads the nemhelper->node_num_map array, node_num_map[i] is the global node number for 00198 // local node number i. 00199 nemhelper->read_node_num_map(); 00200 00201 // The get_cmap_params() function reads in the: 00202 // node_cmap_ids[], 00203 // node_cmap_node_cnts[], 00204 // elem_cmap_ids[], 00205 // elem_cmap_elem_cnts[], 00206 nemhelper->get_cmap_params(); 00207 00208 // Read the IDs of the interior, boundary, and external nodes. This function 00209 // fills the vectors: 00210 // node_mapi[], 00211 // node_mapb[], 00212 // node_mape[] 00213 nemhelper->get_node_map(); 00214 00215 // Read each node communication map for this processor. This function 00216 // fills the vectors of vectors named: 00217 // node_cmap_node_ids[][] 00218 // node_cmap_proc_ids[][] 00219 nemhelper->get_node_cmap(); 00220 00221 libmesh_assert_equal_to (to_uint(nemhelper->num_node_cmaps), nemhelper->node_cmap_node_cnts.size()); 00222 libmesh_assert_equal_to (to_uint(nemhelper->num_node_cmaps), nemhelper->node_cmap_node_ids.size()); 00223 libmesh_assert_equal_to (to_uint(nemhelper->num_node_cmaps), nemhelper->node_cmap_proc_ids.size()); 00224 00225 #ifndef NDEBUG 00226 // We expect the communication maps to be symmetric - e.g. if processor i thinks it 00227 // communicates with processor j, then processor j should also be expecting to 00228 // communicate with i. We can assert that here easily enough with an alltoall, 00229 // but let's only do it when not in optimized mode to limit unnecessary communication. 00230 { 00231 std::vector<unsigned char> pid_send_partner (libMesh::n_processors(), 0); 00232 00233 // strictly speaking, we should expect to communicate with ourself... 00234 pid_send_partner[libMesh::processor_id()] = 1; 00235 00236 // mark each processor id we reference with a node cmap 00237 for (unsigned int cmap=0; cmap<to_uint(nemhelper->num_node_cmaps); cmap++) 00238 { 00239 libmesh_assert_less (to_uint(nemhelper->node_cmap_ids[cmap]), libMesh::n_processors()); 00240 00241 pid_send_partner[nemhelper->node_cmap_ids[cmap]] = 1; 00242 } 00243 00244 // Copy the send pairing so we can catch the receive paring and 00245 // test for equality 00246 const std::vector<unsigned char> pid_recv_partner (pid_send_partner); 00247 00248 CommWorld.alltoall (pid_send_partner); 00249 00250 libmesh_assert (pid_send_partner == pid_recv_partner); 00251 } 00252 #endif 00253 00254 // We now have enough information to infer node ownership. We start by assuming 00255 // we own all the nodes on this processor. We will then interrogate the 00256 // node cmaps and see if a lower-rank processor is associated with any of 00257 // our nodes. If so, then that processor owns the node, not us... 00258 std::vector<unsigned short int> node_ownership (nemhelper->num_internal_nodes + 00259 nemhelper->num_border_nodes, 00260 libMesh::processor_id()); 00261 00262 // a map from processor id to cmap number, to be used later 00263 std::map<unsigned int, unsigned int> pid_to_cmap_map; 00264 00265 // For each node_cmap... 00266 for (unsigned int cmap=0; cmap<to_uint(nemhelper->num_node_cmaps); cmap++) 00267 { 00268 // Good time for error checking... 00269 libmesh_assert_equal_to (to_uint(nemhelper->node_cmap_node_cnts[cmap]), 00270 nemhelper->node_cmap_node_ids[cmap].size()); 00271 00272 libmesh_assert_equal_to (to_uint(nemhelper->node_cmap_node_cnts[cmap]), 00273 nemhelper->node_cmap_proc_ids[cmap].size()); 00274 00275 // In all the samples I have seen, node_cmap_ids[cmap] is the processor 00276 // rank of the remote processor... 00277 const unsigned short int adjcnt_pid_idx = nemhelper->node_cmap_ids[cmap]; 00278 00279 libmesh_assert_less (adjcnt_pid_idx, libMesh::n_processors()); 00280 libmesh_assert_not_equal_to (adjcnt_pid_idx, libMesh::processor_id()); 00281 00282 // We only expect one cmap per adjacent processor 00283 libmesh_assert (!pid_to_cmap_map.count(adjcnt_pid_idx)); 00284 00285 pid_to_cmap_map[adjcnt_pid_idx] = cmap; 00286 00287 // ...and each node in that cmap... 00288 for (unsigned int idx=0; idx<to_uint(nemhelper->node_cmap_node_cnts[cmap]); idx++) 00289 { 00290 // Are the node_cmap_ids and node_cmap_proc_ids really redundant? 00291 libmesh_assert_equal_to (adjcnt_pid_idx, nemhelper->node_cmap_proc_ids[cmap][idx]); 00292 00293 // we are expecting the exodus node numbering to be 1-based... 00294 const unsigned int local_node_idx = nemhelper->node_cmap_node_ids[cmap][idx]-1; 00295 00296 libmesh_assert_less (local_node_idx, node_ownership.size()); 00297 00298 // if the adjacent processor is lower rank than the current 00299 // owner for this node, then it will get the node... 00300 node_ownership[local_node_idx] = 00301 std::min(node_ownership[local_node_idx], adjcnt_pid_idx); 00302 } 00303 } // We now should have established proper node ownership. 00304 00305 // now that ownership is established, we can figure out how many nodes we 00306 // will be responsible for numbering. 00307 unsigned int num_nodes_i_must_number = 0; 00308 00309 for (unsigned int idx=0; idx<node_ownership.size(); idx++) 00310 if (node_ownership[idx] == libMesh::processor_id()) 00311 num_nodes_i_must_number++; 00312 00313 // more error checking... 00314 libmesh_assert_greater_equal (num_nodes_i_must_number, to_uint(nemhelper->num_internal_nodes)); 00315 libmesh_assert (num_nodes_i_must_number <= to_uint(nemhelper->num_internal_nodes + 00316 nemhelper->num_border_nodes)); 00317 if (_verbose) 00318 libMesh::out << "[" << libMesh::processor_id() << "] " 00319 << "num_nodes_i_must_number=" 00320 << num_nodes_i_must_number 00321 << std::endl; 00322 00323 // The call to get_loadbal_param() gets 7 pieces of information. We allgather 00324 // these now across all processors to determine some global numberings. We should 00325 // also gather the number of nodes each processor thinks it will number so that 00326 // we can (i) determine our offset, and (ii) do some error checking. 00327 std::vector<int> all_loadbal_data ( 8 ); 00328 all_loadbal_data[0] = nemhelper->num_internal_nodes; 00329 all_loadbal_data[1] = nemhelper->num_border_nodes; 00330 all_loadbal_data[2] = nemhelper->num_external_nodes; 00331 all_loadbal_data[3] = nemhelper->num_internal_elems; 00332 all_loadbal_data[4] = nemhelper->num_border_elems; 00333 all_loadbal_data[5] = nemhelper->num_node_cmaps; 00334 all_loadbal_data[6] = nemhelper->num_elem_cmaps; 00335 all_loadbal_data[7] = num_nodes_i_must_number; 00336 00337 CommWorld.allgather (all_loadbal_data, /* identical_buffer_sizes = */ true); 00338 00339 // OK, we are now in a position to request new global indices for all the nodes 00340 // we do not own 00341 00342 // Let's get a unique message tag to use for send()/receive() 00343 Parallel::MessageTag nodes_tag = Parallel::Communicator_World.get_unique_tag(12345); 00344 00345 std::vector<std::vector<int> > 00346 needed_node_idxs (nemhelper->num_node_cmaps); // the indices we will ask for 00347 00348 std::vector<Parallel::Request> 00349 needed_nodes_requests (nemhelper->num_node_cmaps); 00350 00351 for (unsigned int cmap=0; cmap<to_uint(nemhelper->num_node_cmaps); cmap++) 00352 { 00353 // We know we will need no more indices than there are nodes 00354 // in this cmap, but that number is an upper bound in general 00355 // since the neighboring processor associated with the cmap 00356 // may not actually own it 00357 needed_node_idxs[cmap].reserve (nemhelper->node_cmap_node_cnts[cmap]); 00358 00359 const unsigned int adjcnt_pid_idx = nemhelper->node_cmap_ids[cmap]; 00360 00361 // ...and each node in that cmap... 00362 for (unsigned int idx=0; idx<to_uint(nemhelper->node_cmap_node_cnts[cmap]); idx++) 00363 { 00364 const unsigned int 00365 local_node_idx = nemhelper->node_cmap_node_ids[cmap][idx]-1, 00366 owning_pid_idx = node_ownership[local_node_idx]; 00367 00368 // add it to the request list for its owning processor. 00369 if (owning_pid_idx == adjcnt_pid_idx) 00370 { 00371 const unsigned int 00372 global_node_idx = nemhelper->node_num_map[local_node_idx]-1; 00373 needed_node_idxs[cmap].push_back(global_node_idx); 00374 } 00375 } 00376 // now post the send for this cmap 00377 CommWorld.send (adjcnt_pid_idx, // destination 00378 needed_node_idxs[cmap], // send buffer 00379 needed_nodes_requests[cmap], // request 00380 nodes_tag); 00381 } // all communication requests for getting updated global indices for border 00382 // nodes have been initiated 00383 00384 // Figure out how many nodes each processor thinks it will number and make sure 00385 // that it adds up to the global number of nodes. Also, set up global node 00386 // index offsets for each processor. 00387 std::vector<unsigned int> 00388 all_num_nodes_i_must_number (libMesh::n_processors()); 00389 00390 for (unsigned int pid=0; pid<libMesh::n_processors(); pid++) 00391 all_num_nodes_i_must_number[pid] = all_loadbal_data[8*pid + 7]; 00392 00393 // The sum of all the entries in this vector should sum to the number of global nodes 00394 libmesh_assert (std::accumulate(all_num_nodes_i_must_number.begin(), 00395 all_num_nodes_i_must_number.end(), 00396 0) == nemhelper->num_nodes_global); 00397 00398 unsigned int my_next_node = 0; 00399 for (unsigned int pid=0; pid<libMesh::processor_id(); pid++) 00400 my_next_node += all_num_nodes_i_must_number[pid]; 00401 00402 const unsigned int my_node_offset = my_next_node; 00403 00404 if (_verbose) 00405 libMesh::out << "[" << libMesh::processor_id() << "] " 00406 << "my_node_offset=" 00407 << my_node_offset 00408 << std::endl; 00409 00410 // Add internal nodes to the ParallelMesh, using the node ID offset we 00411 // computed and the current processor's ID. 00412 for (unsigned int i=0; i<to_uint(nemhelper->num_internal_nodes); ++i) 00413 { 00414 const unsigned int local_node_idx = nemhelper->node_mapi[i]-1; 00415 #ifndef NDEBUG 00416 const unsigned int global_node_idx = nemhelper->node_num_map[local_node_idx]-1; 00417 const unsigned int owning_pid_idx = node_ownership[local_node_idx]; 00418 #endif 00419 00420 // an internal node we do not own? huh?? 00421 libmesh_assert_equal_to (owning_pid_idx, libMesh::processor_id()); 00422 libmesh_assert_less (global_node_idx, to_uint(nemhelper->num_nodes_global)); 00423 00424 // "Catch" the node pointer after addition, make sure the 00425 // ID matches the requested value. 00426 Node* added_node = 00427 mesh.add_point (Point(nemhelper->x[local_node_idx], 00428 nemhelper->y[local_node_idx], 00429 nemhelper->z[local_node_idx]), 00430 my_next_node, 00431 libMesh::processor_id()); 00432 00433 // Make sure the node we added has the ID we thought it would 00434 if (added_node->id() != my_next_node) 00435 { 00436 libMesh::err << "Error, node added with ID " << added_node->id() 00437 << ", but we wanted ID " << my_next_node << std::endl; 00438 } 00439 00440 // update the local->global index map, when we are done 00441 // it will be 0-based. 00442 nemhelper->node_num_map[local_node_idx] = my_next_node++; 00443 } 00444 00445 // Now, for the boundary nodes... We may very well own some of them, 00446 // but there may be others for which we have requested the new global 00447 // id. We expect to be asked for the ids of the ones we own, so 00448 // we need to create a map from the old global id to the new one 00449 // we are about to create. 00450 typedef std::vector<std::pair<unsigned int, unsigned int> > global_idx_mapping_type; 00451 global_idx_mapping_type old_global_to_new_global_map; 00452 old_global_to_new_global_map.reserve (num_nodes_i_must_number // total # i will have 00453 - (my_next_node // amount i have thus far 00454 - my_node_offset)); // this should be exact! 00455 CompareGlobalIdxMappings global_idx_mapping_comp; 00456 00457 for (unsigned int i=0; i<to_uint(nemhelper->num_border_nodes); ++i) 00458 { 00459 const unsigned int 00460 local_node_idx = nemhelper->node_mapb[i]-1, 00461 owning_pid_idx = node_ownership[local_node_idx]; 00462 00463 // if we own it... 00464 if (owning_pid_idx == libMesh::processor_id()) 00465 { 00466 const unsigned int 00467 global_node_idx = nemhelper->node_num_map[local_node_idx]-1; 00468 00469 // we will number it, and create a mapping from its old global index to 00470 // the new global index, for lookup purposes when neighbors come calling 00471 old_global_to_new_global_map.push_back(std::make_pair(global_node_idx, 00472 my_next_node)); 00473 00474 // "Catch" the node pointer after addition, make sure the 00475 // ID matches the requested value. 00476 Node* added_node = 00477 mesh.add_point (Point(nemhelper->x[local_node_idx], 00478 nemhelper->y[local_node_idx], 00479 nemhelper->z[local_node_idx]), 00480 my_next_node, 00481 libMesh::processor_id()); 00482 00483 // Make sure the node we added has the ID we thought it would 00484 if (added_node->id() != my_next_node) 00485 { 00486 libMesh::err << "Error, node added with ID " << added_node->id() 00487 << ", but we wanted ID " << my_next_node << std::endl; 00488 } 00489 00490 // update the local->global index map, when we are done 00491 // it will be 0-based. 00492 nemhelper->node_num_map[local_node_idx] = my_next_node++; 00493 } 00494 } 00495 // That should cover numbering all the nodes which belong to us... 00496 libmesh_assert_equal_to (num_nodes_i_must_number, (my_next_node - my_node_offset)); 00497 00498 // Let's sort the mapping so we can efficiently answer requests 00499 std::sort (old_global_to_new_global_map.begin(), 00500 old_global_to_new_global_map.end(), 00501 global_idx_mapping_comp); 00502 00503 // and it had better be unique... 00504 libmesh_assert (std::unique (old_global_to_new_global_map.begin(), 00505 old_global_to_new_global_map.end(), 00506 global_idx_mapping_equality) == 00507 old_global_to_new_global_map.end()); 00508 00509 // We can now catch incoming requests and process them. for efficiency 00510 // let's do whatever is available next 00511 std::map<unsigned int, std::vector<int> > requested_node_idxs; // the indices asked of us 00512 00513 std::vector<Parallel::Request> requested_nodes_requests(nemhelper->num_node_cmaps); 00514 00515 // We know we will receive the request from a given processor before 00516 // we receive its reply to our request. However, we may receive 00517 // a request and a response from one processor before getting 00518 // a request from another processor. So what we are doing here 00519 // is processing whatever message comes next, while recognizing 00520 // we will receive a request from a processor before receiving 00521 // its reply 00522 std::vector<bool> processed_cmap (nemhelper->num_node_cmaps, false); 00523 00524 for (unsigned int comm_step=0; comm_step<2*to_uint(nemhelper->num_node_cmaps); comm_step++) 00525 { 00526 // query the first message which is available 00527 const Parallel::Status 00528 status (CommWorld.probe (Parallel::any_source, 00529 nodes_tag)); 00530 const unsigned int 00531 requesting_pid_idx = status.source(), 00532 source_pid_idx = status.source(); 00533 00534 // this had better be from a processor we are expecting... 00535 libmesh_assert (pid_to_cmap_map.count(requesting_pid_idx)); 00536 00537 // the local cmap which corresponds to the source processor 00538 const unsigned int cmap = pid_to_cmap_map[source_pid_idx]; 00539 00540 if (!processed_cmap[cmap]) 00541 { 00542 processed_cmap[cmap] = true; 00543 00544 // we should only get one request per paired processor 00545 libmesh_assert (!requested_node_idxs.count(requesting_pid_idx)); 00546 00547 // get a reference to the request buffer for this processor to 00548 // avoid repeated map lookups 00549 std::vector<int> &xfer_buf (requested_node_idxs[requesting_pid_idx]); 00550 00551 // actually receive the message. 00552 CommWorld.receive (requesting_pid_idx, xfer_buf, nodes_tag); 00553 00554 // Fill the request 00555 for (unsigned int i=0; i<xfer_buf.size(); i++) 00556 { 00557 // the requested old global node index, *now 0-based* 00558 const unsigned int old_global_node_idx = xfer_buf[i]; 00559 00560 // find the new global node index for the requested node - 00561 // note that requesting_pid_idx thinks we own this node, 00562 // so we better! 00563 const global_idx_mapping_type::const_iterator it = 00564 std::lower_bound (old_global_to_new_global_map.begin(), 00565 old_global_to_new_global_map.end(), 00566 old_global_node_idx, 00567 global_idx_mapping_comp); 00568 00569 libmesh_assert (it != old_global_to_new_global_map.end()); 00570 libmesh_assert_equal_to (it->first, old_global_node_idx); 00571 libmesh_assert_greater_equal (it->second, my_node_offset); 00572 libmesh_assert_less (it->second, my_next_node); 00573 00574 // overwrite the requested old global node index with the new global index 00575 xfer_buf[i] = it->second; 00576 } 00577 00578 // and send the new global indices back to the processor which asked for them 00579 CommWorld.send (requesting_pid_idx, 00580 xfer_buf, 00581 requested_nodes_requests[cmap], 00582 nodes_tag); 00583 } // done processing the request 00584 00585 // this is the second time we have heard from this processor, 00586 // so it must be its reply to our request 00587 else 00588 { 00589 // a long time ago, we sent off our own requests. now it is time to catch the 00590 // replies and get the new global node numbering. note that for any reply 00591 // we receive, the corresponding nonblocking send from above *must* have been 00592 // completed, since the reply is in response to that request!! 00593 00594 // if we have received a reply, our send *must* have completed 00595 // (note we never actually need to wait on the request) 00596 libmesh_assert (needed_nodes_requests[cmap].test()); 00597 libmesh_assert_equal_to (to_uint(nemhelper->node_cmap_ids[cmap]), source_pid_idx); 00598 00599 // now post the receive for this cmap 00600 CommWorld.receive (source_pid_idx, 00601 needed_node_idxs[cmap], 00602 nodes_tag); 00603 00604 libmesh_assert_less_equal (needed_node_idxs[cmap].size(), 00605 nemhelper->node_cmap_node_ids[cmap].size()); 00606 00607 for (unsigned int i=0,j=0; i<nemhelper->node_cmap_node_ids[cmap].size(); i++) 00608 { 00609 const unsigned int 00610 local_node_idx = nemhelper->node_cmap_node_ids[cmap][i]-1, 00611 owning_pid_idx = node_ownership[local_node_idx]; 00612 00613 // if this node is owned by source_pid_idx, its new global id 00614 // is in the buffer we just received 00615 if (owning_pid_idx == source_pid_idx) 00616 { 00617 libmesh_assert_less (j, needed_node_idxs[cmap].size()); 00618 00619 const unsigned int // now 0-based! 00620 global_node_idx = needed_node_idxs[cmap][j++]; 00621 00622 // "Catch" the node pointer after addition, make sure the 00623 // ID matches the requested value. 00624 Node* added_node = 00625 mesh.add_point (Point(nemhelper->x[local_node_idx], 00626 nemhelper->y[local_node_idx], 00627 nemhelper->z[local_node_idx]), 00628 global_node_idx, 00629 source_pid_idx); 00630 00631 // Make sure the node we added has the ID we thought it would 00632 if (added_node->id() != global_node_idx) 00633 { 00634 libMesh::err << "Error, node added with ID " << added_node->id() 00635 << ", but we wanted ID " << global_node_idx << std::endl; 00636 } 00637 00638 // update the local->global index map, when we are done 00639 // it will be 0-based. 00640 nemhelper->node_num_map[local_node_idx] = global_node_idx; 00641 00642 // we are not really going to use my_next_node again, but we can 00643 // keep incrimenting it to track how many nodes we have added 00644 // to the mesh 00645 my_next_node++; 00646 } 00647 } 00648 } 00649 } // end of node index communication loop 00650 00651 // we had better have added all the nodes we need to! 00652 libmesh_assert_equal_to ((my_next_node - my_node_offset), to_uint(nemhelper->num_nodes)); 00653 00654 // After all that, we should be done with all node-related arrays *except* the 00655 // node_num_map, which we have transformed to use our new numbering... 00656 // So let's clean up the arrays we are done with. 00657 { 00658 Utility::deallocate (nemhelper->node_mapi); 00659 Utility::deallocate (nemhelper->node_mapb); 00660 Utility::deallocate (nemhelper->node_mape); 00661 Utility::deallocate (nemhelper->node_cmap_ids); 00662 Utility::deallocate (nemhelper->node_cmap_node_cnts); 00663 Utility::deallocate (nemhelper->node_cmap_node_ids); 00664 Utility::deallocate (nemhelper->node_cmap_proc_ids); 00665 Utility::deallocate (nemhelper->x); 00666 Utility::deallocate (nemhelper->y); 00667 Utility::deallocate (nemhelper->z); 00668 Utility::deallocate (needed_node_idxs); 00669 Utility::deallocate (node_ownership); 00670 } 00671 00672 Parallel::wait (requested_nodes_requests); 00673 requested_node_idxs.clear(); 00674 00675 // See what the node count is up to now. 00676 if (_verbose) 00677 { 00678 // Report the number of nodes which have been added locally 00679 libMesh::out << "[" << libMesh::processor_id() << "] "; 00680 libMesh::out << "mesh.n_nodes()=" << mesh.n_nodes() << std::endl; 00681 00682 // Reports the number of nodes that have been added in total. 00683 libMesh::out << "[" << libMesh::processor_id() << "] "; 00684 libMesh::out << "mesh.parallel_n_nodes()=" << mesh.parallel_n_nodes() << std::endl; 00685 } 00686 00687 00688 00689 // -------------------------------------------------------------------------------- 00690 // -------------------------------------------------------------------------------- 00691 // -------------------------------------------------------------------------------- 00692 00693 00694 // We can now read in the elements...Exodus stores them in blocks in which all 00695 // elements have the same geometric type. This code is adapted directly from exodusII_io.C 00696 00697 // Assertion: The sum of the border and internal elements on all processors 00698 // should equal nemhelper->num_elems_global 00699 #ifndef NDEBUG 00700 { 00701 int sum_internal_elems=0, sum_border_elems=0; 00702 for (unsigned int j=3,c=0; c<libMesh::n_processors(); j+=8,++c) 00703 sum_internal_elems += all_loadbal_data[j]; 00704 00705 for (unsigned int j=4,c=0; c<libMesh::n_processors(); j+=8,++c) 00706 sum_border_elems += all_loadbal_data[j]; 00707 00708 if (_verbose) 00709 { 00710 libMesh::out << "[" << libMesh::processor_id() << "] "; 00711 libMesh::out << "sum_internal_elems=" << sum_internal_elems << std::endl; 00712 00713 libMesh::out << "[" << libMesh::processor_id() << "] "; 00714 libMesh::out << "sum_border_elems=" << sum_border_elems << std::endl; 00715 } 00716 00717 libmesh_assert_equal_to (sum_internal_elems+sum_border_elems, nemhelper->num_elems_global); 00718 } 00719 #endif 00720 00721 // We need to set the mesh dimension, but the following... 00722 // mesh.set_mesh_dimension(static_cast<unsigned int>(nemhelper->num_dim)); 00723 00724 // ... is not sufficient since some codes report num_dim==3 for two dimensional 00725 // meshes living in 3D, even though all the elements are of 2D type. Therefore, 00726 // we instead use the dimension of the highest element found for the Mesh dimension, 00727 // similar to what is done by the Exodus reader, except here it requires a 00728 // parallel communication. 00729 elems_of_dimension.resize(4, false); // will use 1-based 00730 00731 // Compute my_elem_offset, the amount by which to offset the local elem numbering 00732 // on my processor. 00733 unsigned int my_next_elem = 0; 00734 for (unsigned int pid=0; pid<libMesh::processor_id(); ++pid) 00735 my_next_elem += (all_loadbal_data[8*pid + 3]+ // num_internal_elems, proc pid 00736 all_loadbal_data[8*pid + 4]); // num_border_elems, proc pid 00737 const unsigned int my_elem_offset = my_next_elem; 00738 00739 if (_verbose) 00740 libMesh::out << "[" << libMesh::processor_id() << "] " 00741 << "my_elem_offset=" << my_elem_offset << std::endl; 00742 00743 00744 // Fills in the: 00745 // global_elem_blk_ids[] and 00746 // global_elem_blk_cnts[] arrays. 00747 nemhelper->get_eb_info_global(); 00748 00749 // // Fills in the vectors 00750 // // elem_mapi[num_internal_elems] 00751 // // elem_mapb[num_border_elems ] 00752 // // These tell which of the (locally-numbered) elements are internal and which are border elements. 00753 // // In our test example these arrays are sorted (but non-contiguous), which makes it possible to 00754 // // binary search for each element ID... however I don't think we need to distinguish between the 00755 // // two types, since either can have nodes the boundary! 00756 // nemhelper->get_elem_map(); 00757 00758 // Fills in the vectors of vectors: 00759 // elem_cmap_elem_ids[][] 00760 // elem_cmap_side_ids[][] 00761 // elem_cmap_proc_ids[][] 00762 // These arrays are of size num_elem_cmaps * elem_cmap_elem_cnts[i], i = 0..num_elem_cmaps 00763 nemhelper->get_elem_cmap(); 00764 00765 // Get information about the element blocks: 00766 // (read in the array nemhelper->block_ids[]) 00767 nemhelper->read_block_info(); 00768 00769 // Reads the nemhelper->elem_num_map array, elem_num_map[i] is the global element number for 00770 // local element number i. 00771 nemhelper->read_elem_num_map(); 00772 00773 // Instantiate the ElementMaps interface. This is what translates LibMesh's 00774 // element numbering scheme to Exodus's. 00775 ExodusII_IO_Helper::ElementMaps em; 00776 00777 // Read in the element connectivity for each block by 00778 // looping over all the blocks. 00779 for (unsigned int i=0; i<to_uint(nemhelper->num_elem_blk); i++) 00780 { 00781 // Read the information for block i: For nemhelper->block_ids[i], reads 00782 // elem_type 00783 // num_elem_this_blk 00784 // num_nodes_per_elem 00785 // num_attr 00786 // connect <-- the nodal connectivity array for each element in the block. 00787 nemhelper->read_elem_in_block(i); 00788 00789 // Note that with parallel files it is possible we have no elements in 00790 // this block! 00791 if (!nemhelper->num_elem_this_blk) continue; 00792 00793 // Set subdomain ID based on the block ID. 00794 int subdomain_id = nemhelper->block_ids[i]; 00795 00796 // Create a type string (this uses the null-terminated string ctor). 00797 const std::string type_str ( &(nemhelper->elem_type[0]) ); 00798 00799 // Set any relevant node/edge maps for this element 00800 const ExodusII_IO_Helper::Conversion conv = em.assign_conversion(type_str); 00801 00802 if (_verbose) 00803 libMesh::out << "Reading a block of " << type_str << " elements." << std::endl; 00804 00805 // Loop over all the elements in this block 00806 for (unsigned int j=0; j<to_uint(nemhelper->num_elem_this_blk); j++) 00807 { 00808 Elem* elem = Elem::build (conv.get_canonical_type()).release(); 00809 libmesh_assert (elem); 00810 00811 // Assign subdomain and processor ID to the newly-created Elem. 00812 // Assigning the processor ID beforehand ensures that the Elem is 00813 // not added as an "unpartitioned" element. Note that the element 00814 // numbering in Exodus is also 1-based. 00815 elem->subdomain_id() = subdomain_id; 00816 elem->processor_id() = libMesh::processor_id(); 00817 elem->set_id() = my_next_elem++; 00818 00819 // Mark that we have seen an element of the current element's 00820 // dimension. 00821 elems_of_dimension[elem->dim()] = true; 00822 00823 // Add the created Elem to the Mesh, catch the Elem 00824 // pointer that the Mesh throws back. 00825 elem = mesh.add_elem (elem); 00826 00827 // We are expecting the element "thrown back" by libmesh to have the ID we specified for it... 00828 // Check to see that really is the case. Note that my_next_elem was post-incremented, so 00829 // subtract 1 when performing the check. 00830 if (elem->id() != my_next_elem-1) 00831 { 00832 libMesh::err << "Unexpected ID " 00833 << elem->id() 00834 << " set by parallel mesh. (expecting " 00835 << my_next_elem-1 00836 << ")." << std::endl; 00837 libmesh_error(); 00838 } 00839 00840 // Set all the nodes for this element 00841 if (_verbose) 00842 libMesh::out << "[" << libMesh::processor_id() << "] " 00843 << "Setting nodes for Elem " << elem->id() << std::endl; 00844 00845 for (unsigned int k=0; k<to_uint(nemhelper->num_nodes_per_elem); k++) 00846 { 00847 const unsigned int 00848 gi = (j*nemhelper->num_nodes_per_elem + // index into connectivity array 00849 conv.get_node_map(k)), 00850 local_node_idx = nemhelper->connect[gi]-1, // local node index 00851 global_node_idx = nemhelper->node_num_map[local_node_idx]; // new global node index 00852 00853 // Set node number 00854 elem->set_node(k) = mesh.node_ptr(global_node_idx); 00855 } 00856 } // for (unsigned int j=0; j<nemhelper->num_elem_this_blk; j++) 00857 } // end for (unsigned int i=0; i<nemhelper->num_elem_blk; i++) 00858 00859 libmesh_assert_equal_to ((my_next_elem - my_elem_offset), to_uint(nemhelper->num_elem)); 00860 00861 if (_verbose) 00862 { 00863 // Print local elems_of_dimension information 00864 for (unsigned int i=1; i<elems_of_dimension.size(); ++i) 00865 libMesh::out << "[" << libMesh::processor_id() << "] " 00866 << "elems_of_dimension[" << i << "]=" << elems_of_dimension[i] << std::endl; 00867 } 00868 00869 // Get the max dimension seen on the current processor 00870 unsigned int max_dim_seen = 0; 00871 for (unsigned int i=1; i<elems_of_dimension.size(); ++i) 00872 if (elems_of_dimension[i]) 00873 max_dim_seen = i; 00874 00875 // Do a global max to determine the max dimension seen by all processors. 00876 // It should match -- I don't think we even support calculations on meshes 00877 // with elements of different dimension... 00878 CommWorld.max(max_dim_seen); 00879 00880 if (_verbose) 00881 { 00882 // Print the max element dimension from all processors 00883 libMesh::out << "[" << libMesh::processor_id() << "] " 00884 << "max_dim_seen=" << max_dim_seen << std::endl; 00885 } 00886 00887 // Set the mesh dimension to the largest encountered for an element 00888 mesh.set_mesh_dimension(max_dim_seen); 00889 00890 #if LIBMESH_DIM < 3 00891 if (mesh.mesh_dimension() > LIBMESH_DIM) 00892 { 00893 libMesh::err << "Cannot open dimension " << 00894 mesh.mesh_dimension() << 00895 " mesh file when configured without " << 00896 mesh.mesh_dimension() << "D support." << 00897 std::endl; 00898 libmesh_error(); 00899 } 00900 #endif 00901 00902 00903 // Global sideset information, they are distributed as well, not sure if they will require communication...? 00904 nemhelper->get_ss_param_global(); 00905 00906 if (_verbose) 00907 { 00908 libMesh::out << "[" << libMesh::processor_id() << "] " 00909 << "Read global sideset parameter information." << std::endl; 00910 00911 // These global values should be the same on all processors... 00912 libMesh::out << "[" << libMesh::processor_id() << "] " 00913 << "Number of global sideset IDs: " << nemhelper->global_sideset_ids.size() << std::endl; 00914 } 00915 00916 // Read *local* sideset info the same way it is done in 00917 // exodusII_io_helper. May be called any time after 00918 // nem_helper->read_header(); This sets num_side_sets and resizes 00919 // elem_list, side_list, and id_list to num_elem_all_sidesets. Note 00920 // that there appears to be the same number of sidesets in each file 00921 // but they all have different numbers of entries (some are empty). 00922 // Note that the sum of "nemhelper->num_elem_all_sidesets" over all 00923 // processors should equal the sum of the entries in the "num_global_side_counts" array 00924 // filled up by nemhelper->get_ss_param_global() 00925 nemhelper->read_sideset_info(); 00926 00927 if (_verbose) 00928 { 00929 libMesh::out << "[" << libMesh::processor_id() << "] " 00930 << "nemhelper->num_side_sets = " << nemhelper->num_side_sets << std::endl; 00931 00932 libMesh::out << "[" << libMesh::processor_id() << "] " 00933 << "nemhelper->num_elem_all_sidesets = " << nemhelper->num_elem_all_sidesets << std::endl; 00934 } 00935 00936 #ifdef DEBUG 00937 { 00938 // In DEBUG mode, check that the global number of sidesets reported 00939 // in each nemesis file matches the sum of all local sideset counts 00940 // from each processor. This requires a small communication, so only 00941 // do it in DEBUG mode. 00942 int sum_num_global_side_counts = std::accumulate(nemhelper->num_global_side_counts.begin(), 00943 nemhelper->num_global_side_counts.end(), 00944 0); 00945 00946 // MPI sum up the local files contributions 00947 int sum_num_elem_all_sidesets = nemhelper->num_elem_all_sidesets; 00948 CommWorld.sum(sum_num_elem_all_sidesets); 00949 00950 if (sum_num_global_side_counts != sum_num_elem_all_sidesets) 00951 { 00952 libMesh::err << "Error! global side count reported by Nemesis does not " 00953 << "match the side count reported by the individual files!" << std::endl; 00954 libmesh_error(); 00955 } 00956 } 00957 #endif 00958 00959 // Note that exodus stores sidesets in separate vectors but we want to pack 00960 // them all into a single vector. So when we call read_sideset(), we pass an offset 00961 // into the single vector of all IDs 00962 for (int offset=0, i=0; i<nemhelper->get_num_side_sets(); i++) 00963 { 00964 offset += (i > 0 ? nemhelper->get_num_sides_per_set(i-1) : 0); // Compute new offset 00965 nemhelper->read_sideset (i, offset); 00966 } 00967 00968 // Now that we have the lists of elements, sides, and IDs, we are ready to set them 00969 // in the BoundaryInfo object of our Mesh object. This is slightly different in parallel... 00970 // For example, I think the IDs in each of the split Exodus files are numbered locally, 00971 // and we have to know the appropriate ID for this processor to be able to set the 00972 // entry in BoundaryInfo. This offset should be given by my_elem_offset determined in 00973 // this function... 00974 00975 // Get references to the elem, side, and ID lists for this processor. We should be 00976 // able to get pointers to these elements on the local processor, and add them to the 00977 // BoundaryInfo object 00978 const std::vector<int>& elem_list = nemhelper->get_elem_list(); 00979 const std::vector<int>& side_list = nemhelper->get_side_list(); 00980 const std::vector<int>& id_list = nemhelper->get_id_list(); 00981 00982 00983 // Debugging: 00984 // Print entries of elem_list 00985 // libMesh::out << "[" << libMesh::processor_id() << "] " 00986 // << "elem_list = "; 00987 // for (unsigned int e=0; e<elem_list.size(); e++) 00988 // { 00989 // libMesh::out << elem_list[e] << ", "; 00990 // } 00991 // libMesh::out << std::endl; 00992 00993 // Print entries of side_list 00994 // libMesh::out << "[" << libMesh::processor_id() << "] " 00995 // << "side_list = "; 00996 // for (unsigned int e=0; e<side_list.size(); e++) 00997 // { 00998 // libMesh::out << side_list[e] << ", "; 00999 // } 01000 // libMesh::out << std::endl; 01001 01002 01003 // Loop over the entries of the elem_list, get their pointers from the 01004 // Mesh data structure, and assign the appropriate side to the BoundaryInfo object. 01005 for (unsigned int e=0; e<elem_list.size(); e++) 01006 { 01007 // Calling mesh.elem() feels wrong, for example, in 01008 // ParallelMesh, Mesh::elem() can return NULL so we have to check for this... 01009 // 01010 // Perhaps we should iterate over elements and look them up in 01011 // the elem list instead? Note that the IDs in elem_list are 01012 // not necessarily in order, so if we did instead loop over the 01013 // mesh, we would have to search the (unsorted) elem_list vector 01014 // for each entry! We'll settle for doing some error checking instead. 01015 Elem* elem = mesh.elem(my_elem_offset + (elem_list[e]-1)/*Exodus numbering is 1-based!*/); 01016 01017 if (elem == NULL) 01018 { 01019 libMesh::err << "Mesh returned a NULL pointer when asked for element " 01020 << my_elem_offset << " + " << elem_list[e] << " = " << my_elem_offset+elem_list[e] << std::endl; 01021 libmesh_error(); 01022 } 01023 01024 // The side numberings in libmesh and exodus are not 1:1, so we need to map 01025 // whatever side number is stored in Exodus into a libmesh side number using 01026 // a conv object... 01027 const ExodusII_IO_Helper::Conversion conv = 01028 em.assign_conversion(elem->type()); 01029 01030 // Finally, we are ready to add the element and its side to the BoundaryInfo object. 01031 // Call the version of add_side which takes a pointer, since we have already gone to 01032 // the trouble of getting said pointer... 01033 mesh.boundary_info->add_side (elem, 01034 conv.get_side_map(side_list[e]-1/*Exodus numbering is 1-based*/), 01035 id_list[e]); 01036 } 01037 01038 // Debugging: make sure there are as many boundary conditions in the 01039 // boundary ID object as expected. Note that, at this point, the 01040 // mesh still thinks it's serial, so n_boundary_conds() returns the 01041 // local number of boundary conditions (and is therefore cheap) 01042 // which should match elem_list.size(). 01043 { 01044 std::size_t nbcs = mesh.boundary_info->n_boundary_conds(); 01045 if (nbcs != elem_list.size()) 01046 { 01047 libMesh::err << "[" << libMesh::processor_id() << "] "; 01048 libMesh::err << "BoundaryInfo contains " 01049 << nbcs 01050 << " boundary conditions, while the Exodus file had " 01051 << elem_list.size() 01052 << std::endl; 01053 libmesh_error(); 01054 } 01055 } 01056 01057 // Read global nodeset parameters? We might be able to use this to verify 01058 // something about the local files, but I haven't figured out what yet... 01059 nemhelper->get_ns_param_global(); 01060 01061 // Read local nodeset info 01062 nemhelper->read_nodeset_info(); 01063 01064 if (_verbose) 01065 { 01066 libMesh::out << "[" << libMesh::processor_id() << "] "; 01067 libMesh::out << "nemhelper->num_node_sets=" << nemhelper->num_node_sets << std::endl; 01068 } 01069 01070 // // Debugging, what is currently in nemhelper->node_num_map anyway? 01071 // libMesh::out << "[" << libMesh::processor_id() << "] " 01072 // << "nemhelper->node_num_map = "; 01073 // 01074 // for (unsigned int i=0; i<nemhelper->node_num_map.size(); ++i) 01075 // libMesh::out << nemhelper->node_num_map[i] << ", "; 01076 // libMesh::out << std::endl; 01077 01078 // For each nodeset, 01079 for (int nodeset=0; nodeset<nemhelper->get_num_node_sets(); nodeset++) 01080 { 01081 // Get the user-defined ID associcated with the nodeset 01082 int nodeset_id = nemhelper->get_nodeset_id(nodeset); 01083 01084 if (_verbose) 01085 { 01086 libMesh::out << "[" << libMesh::processor_id() << "] "; 01087 libMesh::out << "nemhelper->get_nodeset_id(" << nodeset << ")=" << nodeset_id << std::endl; 01088 } 01089 01090 // Read the nodeset from file, store them in a vector 01091 nemhelper->read_nodeset(nodeset); 01092 01093 // Get access to that vector 01094 const std::vector<int>& node_list = nemhelper->get_node_list(); 01095 01096 // Add nodes from the node_list to the BoundaryInfo object 01097 for(unsigned int node=0; node<node_list.size(); node++) 01098 { 01099 // Don't run past the end of our node map! 01100 if (to_uint(node_list[node]-1) >= nemhelper->node_num_map.size()) 01101 { 01102 libMesh::err << "Error, index is past the end of node_num_map array!" << std::endl; 01103 libmesh_error(); 01104 } 01105 01106 // We should be able to use the node_num_map data structure set up previously to determine 01107 // the proper global node index. 01108 unsigned global_node_id = nemhelper->node_num_map[ node_list[node]-1 /*Exodus is 1-based!*/ ]; 01109 01110 if (_verbose) 01111 { 01112 libMesh::out << "[" << libMesh::processor_id() << "] " 01113 << "nodeset " << nodeset 01114 << ", local node number: " << node_list[node]-1 01115 << ", global node id: " << global_node_id 01116 << std::endl; 01117 } 01118 01119 // Add the node to the BoundaryInfo object with the proper nodeset_id 01120 mesh.boundary_info->add_node(global_node_id, nodeset_id); 01121 } 01122 } 01123 01124 // See what the elem count is up to now. 01125 if (_verbose) 01126 { 01127 // Report the number of elements which have been added locally 01128 libMesh::out << "[" << libMesh::processor_id() << "] "; 01129 libMesh::out << "mesh.n_elem()=" << mesh.n_elem() << std::endl; 01130 01131 // Reports the number of elements that have been added in total. 01132 libMesh::out << "[" << libMesh::processor_id() << "] "; 01133 libMesh::out << "mesh.parallel_n_elem()=" << mesh.parallel_n_elem() << std::endl; 01134 } 01135 01136 STOP_LOG ("read()","Nemesis_IO"); 01137 01138 // For ParallelMesh, it seems that _is_serial is true by default. A hack to 01139 // make the Mesh think it's parallel might be to call: 01140 mesh.update_post_partitioning(); 01141 mesh.delete_remote_elements(); 01142 01143 // And if that didn't work, then we're actually reading into a 01144 // SerialMesh, so forget about gathering neighboring elements 01145 if (mesh.is_serial()) 01146 return; 01147 01148 // Gather neighboring elements so that the mesh has the proper "ghost" neighbor information. 01149 MeshCommunication().gather_neighboring_elements(libmesh_cast_ref<ParallelMesh&>(mesh)); 01150 }
| void libMesh::MeshInput< MeshBase >::skip_comment_lines | ( | std::istream & | in, | |
| const char | comment_start | |||
| ) | [protected, inherited] |
Reads input from in, skipping all the lines that start with the character comment_start.
Referenced by libMesh::TetGenIO::read(), and libMesh::UCDIO::read_implementation().
| void libMesh::Nemesis_IO::verbose | ( | bool | set_verbosity | ) |
Set the flag indicationg if we should be verbose.
Definition at line 102 of file nemesis_io.C.
References _verbose, nemhelper, and libMesh::ExodusII_IO_Helper::verbose().
| void libMesh::Nemesis_IO::write | ( | const std::string & | base_filename | ) | [virtual] |
This method implements writing a mesh to a specified file.
Implements libMesh::MeshOutput< MeshBase >.
Definition at line 1168 of file nemesis_io.C.
References libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::Nemesis_IO_Helper::create(), libMesh::ExodusII_IO_Helper::ex_err, libMesh::ExodusII_IO_Helper::ex_id, ex_update(), libMesh::Nemesis_IO_Helper::initialize(), libMesh::MeshInput< MeshBase >::mesh(), nemhelper, libMesh::Nemesis_IO_Helper::write_elements(), libMesh::Nemesis_IO_Helper::write_nodal_coordinates(), libMesh::Nemesis_IO_Helper::write_nodesets(), and libMesh::Nemesis_IO_Helper::write_sidesets().
01169 { 01170 // Get a constant reference to the mesh for writing 01171 const MeshBase & mesh = MeshOutput<MeshBase>::mesh(); 01172 01173 // Create the filename for this processor given the base_filename passed in. 01174 std::string nemesis_filename = nemhelper->construct_nemesis_filename(base_filename); 01175 01176 nemhelper->create(nemesis_filename); 01177 01178 // Initialize data structures and write some global Nemesis-specifc data, such as 01179 // communication maps, to file. 01180 nemhelper->initialize(nemesis_filename,mesh); 01181 01182 // Call the Nemesis-specialized version of write_nodal_coordinates() to write 01183 // the nodal coordinates. 01184 nemhelper->write_nodal_coordinates(mesh); 01185 01186 // Call the Nemesis-specialized version of write_elements() to write 01187 // the elements. Note: Must write a zero if a given global block ID has no 01188 // elements... 01189 nemhelper->write_elements(mesh); 01190 01191 // Call our specialized function to write the nodesets 01192 nemhelper->write_nodesets(mesh); 01193 01194 // Call our specialized write_sidesets() function to write the sidesets to file 01195 nemhelper->write_sidesets(mesh); 01196 01197 // Not sure if this is really necessary, but go ahead and flush the file 01198 // once we have written all this stuff. 01199 nemhelper->ex_err = exII::ex_update(nemhelper->ex_id); 01200 }
| virtual void libMesh::MeshOutput< MeshBase >::write_equation_systems | ( | const std::string & | , | |
| const EquationSystems & | , | |||
| const std::set< std::string > * | system_names = NULL | |||
| ) | [virtual, inherited] |
This method implements writing a mesh with data to a specified file where the data is taken from the EquationSystems object.
Referenced by write_timestep().
| void libMesh::Nemesis_IO::write_global_data | ( | const std::vector< Number > & | soln, | |
| const std::vector< std::string > & | names | |||
| ) |
Write out global variables.
Definition at line 1292 of file nemesis_io.C.
References _timestep, libMesh::ExodusII_IO_Helper::created(), libMesh::err, libMesh::ExodusII_IO_Helper::initialize_global_variables(), nemhelper, and libMesh::ExodusII_IO_Helper::write_global_values().
01294 { 01295 if (!nemhelper->created()) 01296 { 01297 libMesh::err << "ERROR, Nemesis file must be initialized " 01298 << "before outputting global variables.\n" 01299 << std::endl; 01300 libmesh_error(); 01301 } 01302 01303 // Call the Exodus writer implementation 01304 nemhelper->initialize_global_variables( names ); 01305 nemhelper->write_global_values( soln, _timestep); 01306 }
| void libMesh::Nemesis_IO::write_information_records | ( | const std::vector< std::string > & | records | ) |
Write out information records.
Definition at line 1324 of file nemesis_io.C.
References libMesh::ExodusII_IO_Helper::created(), libMesh::err, nemhelper, and libMesh::ExodusII_IO_Helper::write_information_records().
01325 { 01326 if (!nemhelper->created()) 01327 { 01328 libMesh::err << "ERROR, Nemesis file must be initialized " 01329 << "before outputting information records.\n" 01330 << std::endl; 01331 libmesh_error(); 01332 } 01333 01334 // Call the Exodus writer implementation 01335 nemhelper->write_information_records( records ); 01336 }
| void libMesh::Nemesis_IO::write_nodal_data | ( | const std::string & | fname, | |
| const std::vector< Number > & | soln, | |||
| const std::vector< std::string > & | names | |||
| ) | [virtual] |
Output a nodal solution.
Reimplemented from libMesh::MeshOutput< MeshBase >.
Definition at line 1241 of file nemesis_io.C.
References _timestep, libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::Nemesis_IO_Helper::create(), libMesh::ExodusII_IO_Helper::created(), libMesh::Nemesis_IO_Helper::initialize(), libMesh::ExodusII_IO_Helper::initialize_nodal_variables(), libMesh::MeshInput< MeshBase >::mesh(), nemhelper, libMesh::Nemesis_IO_Helper::write_elements(), libMesh::Nemesis_IO_Helper::write_nodal_coordinates(), libMesh::Nemesis_IO_Helper::write_nodal_solution(), libMesh::Nemesis_IO_Helper::write_nodesets(), and libMesh::Nemesis_IO_Helper::write_sidesets().
01244 { 01245 START_LOG("write_nodal_data()", "Nemesis_IO"); 01246 01247 const MeshBase & mesh = MeshOutput<MeshBase>::mesh(); 01248 01249 std::string nemesis_filename = nemhelper->construct_nemesis_filename(base_filename); 01250 01251 if (!nemhelper->created()) 01252 { 01253 nemhelper->create(nemesis_filename); 01254 nemhelper->initialize(nemesis_filename,mesh); 01255 nemhelper->write_nodal_coordinates(mesh); 01256 nemhelper->write_elements(mesh); 01257 nemhelper->write_nodesets(mesh); 01258 nemhelper->write_sidesets(mesh); 01259 01260 // If we don't have any nodes written out on this processor, 01261 // Exodus seems to like us better if we don't try to write out any 01262 // variable names too... 01263 nemhelper->initialize_nodal_variables(names); 01264 } 01265 01266 nemhelper->write_nodal_solution(soln, names, _timestep); 01267 01268 STOP_LOG("write_nodal_data()", "Nemesis_IO"); 01269 }
| void libMesh::Nemesis_IO::write_timestep | ( | const std::string & | fname, | |
| const EquationSystems & | es, | |||
| const int | timestep, | |||
| const Real | time | |||
| ) |
Write one timestep's worth of the solution.
Definition at line 1215 of file nemesis_io.C.
References _timestep, nemhelper, libMesh::MeshOutput< MeshBase >::write_equation_systems(), and libMesh::ExodusII_IO_Helper::write_timestep().
01219 { 01220 _timestep=timestep; 01221 write_equation_systems(fname,es); 01222 01223 nemhelper->write_timestep(timestep, time); 01224 }
Member Data Documentation
const bool libMesh::MeshOutput< MeshBase >::_is_parallel_format [protected, inherited] |
Flag specifying whether this format is parallel-capable. If this is false (default) I/O is only permitted when the mesh has been serialized.
Definition at line 126 of file mesh_output.h.
Referenced by libMesh::PostscriptIO::write(), libMesh::FroIO::write(), libMesh::EnsightIO::write(), and libMesh::DivaIO::write().
int libMesh::Nemesis_IO::_timestep [private] |
Definition at line 112 of file nemesis_io.h.
Referenced by write_global_data(), write_nodal_data(), and write_timestep().
bool libMesh::Nemesis_IO::_verbose [private] |
Definition at line 114 of file nemesis_io.h.
std::vector<bool> libMesh::MeshInput< MeshBase >::elems_of_dimension [protected, inherited] |
A vector of bools describing what dimension elements have been encountered when reading a mesh.
Definition at line 93 of file mesh_input.h.
Referenced by libMesh::GMVIO::_read_one_cell(), libMesh::UNVIO::element_in(), libMesh::VTKIO::read(), read(), libMesh::GMVIO::read(), libMesh::ExodusII_IO::read(), libMesh::UNVIO::read_implementation(), libMesh::UCDIO::read_implementation(), libMesh::LegacyXdrIO::read_mesh(), and libMesh::XdrIO::read_serialized_connectivity().
Nemesis_IO_Helper* libMesh::Nemesis_IO::nemhelper [private] |
Definition at line 110 of file nemesis_io.h.
Referenced by read(), verbose(), write(), write_global_data(), write_information_records(), write_nodal_data(), write_timestep(), and ~Nemesis_IO().
The documentation for this class was generated from the following files:
Site Created By: libMesh Developers
Last modified: February 05 2013 19:55:30 UTC
Hosted By: