libMesh::GMVIO Class Reference

#include <gmv_io.h>

Inheritance diagram for libMesh::GMVIO:

Public Member Functions

 GMVIO (const MeshBase &)
 
 GMVIO (MeshBase &)
 
virtual void write (const std::string &) libmesh_override
 
virtual void read (const std::string &mesh_file) libmesh_override
 
virtual void write_nodal_data (const std::string &, const std::vector< Number > &, const std::vector< std::string > &) libmesh_override
 
bool & binary ()
 
bool & discontinuous ()
 
bool & partitioning ()
 
bool & write_subdomain_id_as_material ()
 
bool & subdivide_second_order ()
 
bool & p_levels ()
 
void write_discontinuous_gmv (const std::string &name, const EquationSystems &es, const bool write_partitioning, const std::set< std::string > *system_names=libmesh_nullptr) const
 
void write_ascii_new_impl (const std::string &, const std::vector< Number > *=libmesh_nullptr, const std::vector< std::string > *=libmesh_nullptr)
 
void add_cell_centered_data (const std::string &cell_centered_data_name, const std::vector< Real > *cell_centered_data_vals)
 
void copy_nodal_solution (EquationSystems &es)
 
virtual void write_equation_systems (const std::string &, const EquationSystems &, const std::set< std::string > *system_names=libmesh_nullptr)
 
virtual void write_nodal_data (const std::string &, const NumericVector< Number > &, const std::vector< std::string > &)
 
unsigned int & ascii_precision ()
 

Protected Member Functions

MeshBasemesh ()
 
void set_n_partitions (unsigned int n_parts)
 
void skip_comment_lines (std::istream &in, const char comment_start)
 
const MeshBasemesh () const
 

Protected Attributes

std::vector< bool > elems_of_dimension
 
const bool _is_parallel_format
 
const bool _serial_only_needed_on_proc_0
 

Private Member Functions

void write_ascii_old_impl (const std::string &, const std::vector< Number > *=libmesh_nullptr, const std::vector< std::string > *=libmesh_nullptr)
 
void write_binary (const std::string &, const std::vector< Number > *=libmesh_nullptr, const std::vector< std::string > *=libmesh_nullptr)
 
void _read_nodes ()
 
void _read_one_cell ()
 
ElemType gmv_elem_to_libmesh_elem (std::string elemname)
 
void _read_materials ()
 
void _read_var ()
 

Static Private Member Functions

static std::map< std::string, ElemTypebuild_reading_element_map ()
 

Private Attributes

bool _binary
 
bool _discontinuous
 
bool _partitioning
 
bool _write_subdomain_id_as_material
 
bool _subdivide_second_order
 
bool _p_levels
 
std::map< std::string, const std::vector< Real > * > _cell_centered_data
 
unsigned int _next_elem_id
 
std::map< std::string, std::vector< Number > > _nodal_data
 

Static Private Attributes

static std::map< std::string, ElemType_reading_element_map = GMVIO::build_reading_element_map()
 

Detailed Description

This class implements writing meshes in the GMV format. For a full description of the GMV format and to obtain the GMV software see http://www.generalmeshviewer.com

Author
Benjamin S. Kirk
Date
2004

Definition at line 46 of file gmv_io.h.

Constructor & Destructor Documentation

libMesh::GMVIO::GMVIO ( const MeshBase mesh)
explicit

Constructor. Takes a reference to a constant mesh object. This constructor will only allow us to write the mesh.

Definition at line 238 of file gmv_io.C.

238  :
240  _binary (false),
241  _discontinuous (false),
242  _partitioning (true),
245  _p_levels (true),
246  _next_elem_id (0)
247 {
248 }
bool _binary
Definition: gmv_io.h:191
unsigned int _next_elem_id
Definition: gmv_io.h:231
bool _discontinuous
Definition: gmv_io.h:196
bool _write_subdomain_id_as_material
Definition: gmv_io.h:207
bool _subdivide_second_order
Definition: gmv_io.h:212
bool _p_levels
Definition: gmv_io.h:217
bool _partitioning
Definition: gmv_io.h:201
libMesh::GMVIO::GMVIO ( MeshBase mesh)
explicit

Constructor. Takes a writable reference to a mesh object. This constructor is required to let us read in a mesh.

Definition at line 252 of file gmv_io.C.

252  :
254  MeshOutput<MeshBase>(mesh),
255  _binary (false),
256  _discontinuous (false),
257  _partitioning (true),
260  _p_levels (true),
261  _next_elem_id (0)
262 {
263 }
bool _binary
Definition: gmv_io.h:191
unsigned int _next_elem_id
Definition: gmv_io.h:231
bool _discontinuous
Definition: gmv_io.h:196
bool _write_subdomain_id_as_material
Definition: gmv_io.h:207
bool _subdivide_second_order
Definition: gmv_io.h:212
bool _p_levels
Definition: gmv_io.h:217
bool _partitioning
Definition: gmv_io.h:201

Member Function Documentation

void libMesh::GMVIO::_read_materials ( )
private

Definition at line 2191 of file gmv_io.C.

References libMesh::MeshInput< MT >::mesh().

Referenced by read().

2192 {
2193 #ifdef LIBMESH_HAVE_GMV
2194 
2195  // LibMesh assigns materials on a per-cell basis
2196  libmesh_assert_equal_to (GMVLib::gmv_data.datatype, CELL);
2197 
2198  // // Material names: LibMesh has no use for these currently...
2199  // libMesh::out << "Number of material names="
2200  // << GMVLib::gmv_data.num
2201  // << std::endl;
2202 
2203  // for (int i = 0; i < GMVLib::gmv_data.num; i++)
2204  // {
2205  // // Build a 32-char string from the appropriate entries
2206  // std::string mat_string(&GMVLib::gmv_data.chardata1[i*33], 32);
2207 
2208  // libMesh::out << "Material name " << i << ": " << mat_string << std::endl;
2209  // }
2210 
2211  // // Material labels: These correspond to (1-based) CPU IDs, and
2212  // // there should be 1 of these for each element.
2213  // libMesh::out << "Number of material labels = "
2214  // << GMVLib::gmv_data.nlongdata1
2215  // << std::endl;
2216 
2217  for (int i = 0; i < GMVLib::gmv_data.nlongdata1; i++)
2218  {
2219  // Debugging Info
2220  // libMesh::out << "Material ID " << i << ": "
2221  // << GMVLib::gmv_data.longdata1[i]
2222  // << std::endl;
2223 
2224  MeshInput<MeshBase>::mesh().elem_ref(i).processor_id() =
2225  cast_int<processor_id_type>(GMVLib::gmv_data.longdata1[i]-1);
2226  }
2227 
2228 #endif
2229 }
void libMesh::GMVIO::_read_nodes ( )
private

Helper functions for reading nodes/cells from a GMV file

Definition at line 2234 of file gmv_io.C.

References libMesh::MeshInput< MT >::mesh().

Referenced by read().

2235 {
2236 #ifdef LIBMESH_HAVE_GMV
2237  // Debugging
2238  // libMesh::out << "gmv_data.datatype = " << GMVLib::gmv_data.datatype << std::endl;
2239 
2240  // LibMesh writes UNSTRUCT=100 node data
2241  libmesh_assert_equal_to (GMVLib::gmv_data.datatype, UNSTRUCT);
2242 
2243  // The nodal data is stored in gmv_data.doubledata{1,2,3}
2244  // and is nnodes long
2245  for (int i = 0; i < GMVLib::gmv_data.num; i++)
2246  {
2247  // Debugging
2248  // libMesh::out << "(x,y,z)="
2249  // << "("
2250  // << GMVLib::gmv_data.doubledata1[i] << ","
2251  // << GMVLib::gmv_data.doubledata2[i] << ","
2252  // << GMVLib::gmv_data.doubledata3[i]
2253  // << ")"
2254  // << std::endl;
2255 
2256  // Add the point to the Mesh
2257  MeshInput<MeshBase>::mesh().add_point(Point(GMVLib::gmv_data.doubledata1[i],
2258  GMVLib::gmv_data.doubledata2[i],
2259  GMVLib::gmv_data.doubledata3[i]), i);
2260  }
2261 #endif
2262 }
A geometric point in (x,y,z) space.
Definition: point.h:38
void libMesh::GMVIO::_read_one_cell ( )
private

Definition at line 2265 of file gmv_io.C.

References _next_elem_id, libMesh::MeshBase::add_elem(), libMesh::Elem::build(), libMesh::Elem::dim(), libMesh::MeshInput< MeshBase >::elems_of_dimension, gmv_elem_to_libmesh_elem(), libMesh::libmesh_assert(), libMesh::MeshInput< MeshBase >::mesh(), libMesh::MeshInput< MT >::mesh(), libMesh::MeshBase::node_ptr(), libMesh::DofObject::set_id(), and libMesh::Elem::set_node().

Referenced by read().

2266 {
2267 #ifdef LIBMESH_HAVE_GMV
2268  // Debugging
2269  // libMesh::out << "gmv_data.datatype=" << GMVLib::gmv_data.datatype << std::endl;
2270 
2271  // This is either a REGULAR=111 cell or
2272  // the ENDKEYWORD=207 of the cells
2273 #ifndef NDEBUG
2274  bool recognized =
2275  (GMVLib::gmv_data.datatype==REGULAR) ||
2276  (GMVLib::gmv_data.datatype==ENDKEYWORD);
2277 #endif
2278  libmesh_assert (recognized);
2279 
2281 
2282  if (GMVLib::gmv_data.datatype == REGULAR)
2283  {
2284  // Debugging
2285  // libMesh::out << "Name of the cell is: " << GMVLib::gmv_data.name1 << std::endl;
2286  // libMesh::out << "Cell has " << GMVLib::gmv_data.num2 << " vertices." << std::endl;
2287 
2288  // We need a mapping from GMV element types to LibMesh
2289  // ElemTypes. Basically the reverse of the eletypes
2290  // std::map above.
2291  //
2292  // FIXME: Since Quad9's apparently don't exist for GMV, and since
2293  // In general we write linear sub-elements to GMV files, we need
2294  // to be careful to read back in exactly what we wrote out...
2295  //
2296  // The gmv_data.name1 field is padded with blank characters when
2297  // it is read in by GMV, so the function that converts it to the
2298  // libmesh element type needs to take that into account.
2299  ElemType type = this->gmv_elem_to_libmesh_elem(GMVLib::gmv_data.name1);
2300 
2301  Elem * elem = Elem::build(type).release();
2302  elem->set_id(_next_elem_id++);
2303 
2304  // Get the ElementDefinition object for this element type
2305  const ElementDefinition & eledef = eletypes[type];
2306 
2307  // Print out the connectivity information for
2308  // this cell.
2309  for (int i=0; i<GMVLib::gmv_data.num2; i++)
2310  {
2311  // Debugging
2312  // libMesh::out << "Vertex " << i << " is node " << GMVLib::gmv_data.longdata1[i] << std::endl;
2313 
2314  // Map index i to GMV's numbering scheme
2315  unsigned mapped_i = eledef.node_map[i];
2316 
2317  // Note: Node numbers (as stored in libmesh) are 1-based
2318  elem->set_node(i) = mesh.node_ptr
2319  (cast_int<dof_id_type>(GMVLib::gmv_data.longdata1[mapped_i]-1));
2320  }
2321 
2322  elems_of_dimension[elem->dim()] = true;
2323 
2324  // Add the newly-created element to the mesh
2325  mesh.add_elem(elem);
2326  }
2327 
2328 
2329  if (GMVLib::gmv_data.datatype == ENDKEYWORD)
2330  {
2331  // There isn't a cell to read, so we just return
2332  return;
2333  }
2334 
2335 #endif
2336 }
virtual Node *& set_node(const unsigned int i)
Definition: elem.h:1789
static UniquePtr< Elem > build(const ElemType type, Elem *p=libmesh_nullptr)
Definition: elem.C:234
std::vector< bool > elems_of_dimension
Definition: mesh_input.h:97
unsigned int _next_elem_id
Definition: gmv_io.h:231
The base class for all geometric element types.
Definition: elem.h:86
virtual const Node * node_ptr(const dof_id_type i) const =0
Base class for Mesh.
Definition: mesh_base.h:67
dof_id_type & set_id()
Definition: dof_object.h:633
libmesh_assert(j)
ElemType gmv_elem_to_libmesh_elem(std::string elemname)
Definition: gmv_io.C:2339
virtual Elem * add_elem(Elem *e)=0
virtual unsigned int dim() const =0
void libMesh::GMVIO::_read_var ( )
private

Definition at line 2179 of file gmv_io.C.

References _nodal_data.

Referenced by read().

2180 {
2181 #ifdef LIBMESH_HAVE_GMV
2182 
2183  // Copy all the variable's values into a local storage vector.
2184  _nodal_data.insert ( std::make_pair(std::string(GMVLib::gmv_data.name1),
2185  std::vector<Number>(GMVLib::gmv_data.doubledata1, GMVLib::gmv_data.doubledata1+GMVLib::gmv_data.num) ) );
2186 #endif
2187 }
std::map< std::string, std::vector< Number > > _nodal_data
Definition: gmv_io.h:236
void libMesh::GMVIO::add_cell_centered_data ( const std::string &  cell_centered_data_name,
const std::vector< Real > *  cell_centered_data_vals 
)

Takes a vector of cell-centered data to be plotted. You must ensure that for every active element e, v[e->id()] is a valid number. You can add an arbitrary number of different cell-centered data sets by calling this function multiple times.

.) GMV does not like spaces in the cell_centered_data_name .) No matter what order you add cell-centered data, it will be output alphabetically.

Definition at line 2005 of file gmv_io.C.

References _cell_centered_data, and libMesh::libmesh_assert().

Referenced by p_levels().

2007 {
2008  libmesh_assert(cell_centered_data_vals);
2009 
2010  // Make sure there are *at least* enough entries for all the active elements.
2011  // There can also be entries for inactive elements, they will be ignored.
2012  // libmesh_assert_greater_equal (cell_centered_data_vals->size(),
2013  // MeshOutput<MeshBase>::mesh().n_active_elem());
2014  this->_cell_centered_data[cell_centered_data_name] = cell_centered_data_vals;
2015 }
libmesh_assert(j)
std::map< std::string, const std::vector< Real > * > _cell_centered_data
Definition: gmv_io.h:225
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(), write_ascii_new_impl(), and write_ascii_old_impl().

bool& libMesh::GMVIO::binary ( )
inline

Flag indicating whether or not to write a binary file. While binary files may end up being smaller than equivalent ASCII files, they are harder to debug if anything goes wrong, since they are not human-readable.

Definition at line 95 of file gmv_io.h.

References _binary.

Referenced by write(), and write_nodal_data().

95 { return _binary; }
bool _binary
Definition: gmv_io.h:191
std::map< std::string, ElemType > libMesh::GMVIO::build_reading_element_map ( )
staticprivate

Static function used to build the _reading_element_map.

Definition at line 205 of file gmv_io.C.

References libMesh::EDGE2, libMesh::EDGE3, libMesh::HEX20, libMesh::HEX27, libMesh::HEX8, libMesh::PRISM15, libMesh::PRISM6, libMesh::QUAD4, libMesh::QUAD8, libMesh::TET10, libMesh::TET4, libMesh::TRI3, and libMesh::TRI6.

206 {
207  std::map<std::string, ElemType> ret;
208 
209  ret["line"] = EDGE2;
210  ret["tri"] = TRI3;
211  ret["tri3"] = TRI3;
212  ret["quad"] = QUAD4;
213  ret["tet"] = TET4;
214  ret["ptet4"] = TET4;
215  ret["hex"] = HEX8;
216  ret["phex8"] = HEX8;
217  ret["prism"] = PRISM6;
218  ret["pprism6"] = PRISM6;
219  ret["phex20"] = HEX20;
220  ret["phex27"] = HEX27;
221  ret["pprism15"] = PRISM15;
222  ret["ptet10"] = TET10;
223  ret["6tri"] = TRI6;
224  ret["8quad"] = QUAD8;
225  ret["3line"] = EDGE3;
226 
227  // Unsupported/Unused types
228  // ret["vface2d"]
229  // ret["vface3d"]
230  // ret["pyramid"]
231  // ret["ppyrmd5"]
232  // ret["ppyrmd13"]
233 
234  return ret;
235 }
void libMesh::GMVIO::copy_nodal_solution ( EquationSystems es)

If we read in a nodal solution while reading in a mesh, we can attempt to copy that nodal solution into an EquationSystems object.

Definition at line 2356 of file gmv_io.C.

References _nodal_data, end, libMesh::err, libMesh::FEType::family, libMesh::FIRST, libMesh::OrderWrapper::get_order(), libMesh::EquationSystems::get_system(), libMesh::System::has_variable(), libMesh::LAGRANGE, libMesh::libmesh_assert(), libMesh::MeshInput< MT >::mesh(), libMesh::EquationSystems::n_systems(), libMesh::FEType::order, libMesh::System::solution, libMesh::sys, libMesh::System::update(), libMesh::System::variable_number(), and libMesh::System::variable_type().

Referenced by p_levels().

2357 {
2358  // Check for easy return if there isn't any nodal data
2359  if (_nodal_data.empty())
2360  {
2361  libMesh::err << "Unable to copy nodal solution: No nodal "
2362  << "solution has been read in from file." << std::endl;
2363  return;
2364  }
2365 
2366  // Be sure there is at least one system
2367  libmesh_assert (es.n_systems());
2368 
2369  // Keep track of variable names which have been found and
2370  // copied already. This could be used to prevent us from
2371  // e.g. copying the same var into 2 different systems ...
2372  // but this seems unlikely. Also, it is used to tell if
2373  // any variables which were read in were not successfully
2374  // copied to the EquationSystems.
2375  std::set<std::string> vars_copied;
2376 
2377  // For each entry in the nodal data map, try to find a system
2378  // that has the same variable key name.
2379  for (unsigned int sys=0; sys<es.n_systems(); ++sys)
2380  {
2381  // Get a generic refernence to the current System
2382  System & system = es.get_system(sys);
2383 
2384  // And a reference to that system's dof_map
2385  // const DofMap & dof_map = system.get_dof_map();
2386 
2387  // For each var entry in the _nodal_data map, try to find
2388  // that var in the system
2389  std::map<std::string, std::vector<Number> >::iterator it = _nodal_data.begin();
2390  const std::map<std::string, std::vector<Number> >::iterator end = _nodal_data.end();
2391  for (; it != end; ++it)
2392  {
2393  std::string var_name = it->first;
2394  // libMesh::out << "Searching for var " << var_name << " in system " << sys << std::endl;
2395 
2396  if (system.has_variable(var_name))
2397  {
2398  // Check if there are as many nodes in the mesh as there are entries
2399  // in the stored nodal data vector
2400  libmesh_assert_equal_to ( it->second.size(), MeshInput<MeshBase>::mesh().n_nodes() );
2401 
2402  const unsigned int var_num = system.variable_number(var_name);
2403 
2404  // libMesh::out << "Variable "
2405  // << var_name
2406  // << " is variable "
2407  // << var_num
2408  // << " in system " << sys << std::endl;
2409 
2410  // The only type of nodal data we can read in from GMV is for
2411  // linear LAGRANGE type elements.
2412  const FEType & fe_type = system.variable_type(var_num);
2413  if ((fe_type.order.get_order() != FIRST) || (fe_type.family != LAGRANGE))
2414  {
2415  libMesh::err << "Only FIRST-order LAGRANGE variables can be read from GMV files. "
2416  << "Skipping variable " << var_name << std::endl;
2417  break;
2418  }
2419 
2420 
2421  // Loop over the stored vector's entries, inserting them into
2422  // the System's solution if appropriate.
2423  for (std::size_t i=0; i<it->second.size(); ++i)
2424  {
2425  // Since this var came from a GMV file, the index i corresponds to
2426  // the (single) DOF value of the current variable for node i.
2427  const unsigned int dof_index =
2428  MeshInput<MeshBase>::mesh().node_ptr(i)->dof_number(sys, // system #
2429  var_num, // var #
2430  0); // component #, always zero for LAGRANGE
2431 
2432  // libMesh::out << "Value " << i << ": "
2433  // << it->second [i]
2434  // << ", dof index="
2435  // << dof_index << std::endl;
2436 
2437  // If the dof_index is local to this processor, set the value
2438  if ((dof_index >= system.solution->first_local_index()) &&
2439  (dof_index < system.solution->last_local_index()))
2440  system.solution->set (dof_index, it->second [i]);
2441  } // end loop over my GMVIO's copy of the solution
2442 
2443  // Add the most recently copied var to the set of copied vars
2444  vars_copied.insert (var_name);
2445  } // end if (system.has_variable)
2446  } // end for loop over _nodal_data
2447 
2448  // Communicate parallel values before going to the next system.
2449  system.solution->close();
2450  system.update();
2451 
2452  } // end loop over all systems
2453 
2454 
2455 
2456  // Warn the user if any GMV variables were not successfully copied over to the EquationSystems object
2457  {
2458  std::map<std::string, std::vector<Number> >::iterator it = _nodal_data.begin();
2459  const std::map<std::string, std::vector<Number> >::iterator end = _nodal_data.end();
2460 
2461  for (; it != end; ++it)
2462  {
2463  if (vars_copied.find( it->first ) == vars_copied.end())
2464  {
2465  libMesh::err << "Warning: Variable "
2466  << it->first
2467  << " was not copied to the EquationSystems object."
2468  << std::endl;
2469  }
2470  }
2471  }
2472 
2473 }
Manages the family, order, etc. parameters for a given FE.
Definition: fe_type.h:178
FEFamily family
Definition: fe_type.h:203
int get_order() const
Definition: fe_type.h:77
ImplicitSystem & sys
bool has_variable(const std::string &var) const
Definition: system.C:1250
OrderWrapper order
Definition: fe_type.h:197
IterBase * end
const FEType & variable_type(const unsigned int i) const
Definition: system.h:2162
std::map< std::string, std::vector< Number > > _nodal_data
Definition: gmv_io.h:236
libmesh_assert(j)
unsigned short int variable_number(const std::string &var) const
Definition: system.C:1257
UniquePtr< NumericVector< Number > > solution
Definition: system.h:1521
OStreamProxy err(std::cerr)
unsigned int n_systems() const
virtual void update()
Definition: system.C:420
const T_sys & get_system(const std::string &name) const
bool& libMesh::GMVIO::discontinuous ( )
inline

Flag indicating whether or not to write the mesh as discontinuous cell patches

Definition at line 101 of file gmv_io.h.

References _discontinuous.

101 { return _discontinuous; }
bool _discontinuous
Definition: gmv_io.h:196
ElemType libMesh::GMVIO::gmv_elem_to_libmesh_elem ( std::string  elemname)
private

Definition at line 2339 of file gmv_io.C.

References _reading_element_map.

Referenced by _read_one_cell().

2340 {
2341  // Erase any whitespace padding in name coming from gmv before performing comparison.
2342  elemname.erase(std::remove_if(elemname.begin(), elemname.end(), isspace), elemname.end());
2343 
2344  // Look up the string in our string->ElemType name.
2345  std::map<std::string, ElemType>::iterator it = _reading_element_map.find(elemname);
2346 
2347  if (it == _reading_element_map.end())
2348  libmesh_error_msg("Uknown/unsupported element: " << elemname << " was read.");
2349 
2350  return it->second;
2351 }
static std::map< std::string, ElemType > _reading_element_map
Definition: gmv_io.h:242
MeshBase & libMesh::MeshInput< MeshBase >::mesh ( )
protectedinherited
Returns
The object as a writable reference.

Referenced by _read_one_cell(), libMesh::TetGenIO::element_in(), libMesh::UNVIO::elements_in(), libMesh::UNVIO::elements_out(), libMesh::UNVIO::groups_in(), libMesh::TetGenIO::node_in(), libMesh::UNVIO::nodes_in(), libMesh::UNVIO::nodes_out(), read(), libMesh::ExodusII_IO::read(), libMesh::CheckpointIO::read(), libMesh::XdrIO::read(), libMesh::CheckpointIO::read_bcs(), libMesh::CheckpointIO::read_connectivity(), libMesh::UCDIO::read_implementation(), libMesh::UNVIO::read_implementation(), libMesh::GmshIO::read_mesh(), libMesh::CheckpointIO::read_nodes(), libMesh::CheckpointIO::read_nodesets(), libMesh::CheckpointIO::read_remote_elem(), 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::OFFIO::read_stream(), libMesh::MatlabIO::read_stream(), libMesh::CheckpointIO::read_subdomain_names(), libMesh::UCDIO::UCDIO(), libMesh::VTKIO::VTKIO(), libMesh::TetGenIO::write(), libMesh::ExodusII_IO::write(), libMesh::CheckpointIO::write(), libMesh::XdrIO::write(), write_ascii_new_impl(), write_ascii_old_impl(), libMesh::CheckpointIO::write_bcs(), write_binary(), write_discontinuous_gmv(), libMesh::ExodusII_IO::write_element_data(), libMesh::UCDIO::write_implementation(), libMesh::GmshIO::write_mesh(), libMesh::UCDIO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data_common(), libMesh::ExodusII_IO::write_nodal_data_discontinuous(), libMesh::CheckpointIO::write_nodesets(), libMesh::XdrIO::write_parallel(), libMesh::GmshIO::write_post(), 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(), and libMesh::CheckpointIO::write_subdomain_names().

bool& libMesh::GMVIO::p_levels ( )
inline

Flag indicating whether or not to write p level information for p refined meshes

Definition at line 126 of file gmv_io.h.

References _p_levels, add_cell_centered_data(), copy_nodal_solution(), libmesh_nullptr, libMesh::Quality::name(), write_ascii_new_impl(), write_ascii_old_impl(), write_binary(), and write_discontinuous_gmv().

Referenced by write_ascii_new_impl(), write_ascii_old_impl(), and write_binary().

126 { return _p_levels; }
bool _p_levels
Definition: gmv_io.h:217
bool& libMesh::GMVIO::partitioning ( )
inline

Flag indicating whether or not to write the partitioning information for the mesh.

Definition at line 107 of file gmv_io.h.

References _partitioning.

Referenced by libMesh::NameBasedIO::write(), write_ascii_new_impl(), write_ascii_old_impl(), write_binary(), and libMesh::NameBasedIO::write_nodal_data().

107 { return _partitioning; }
bool _partitioning
Definition: gmv_io.h:201
void libMesh::GMVIO::read ( const std::string &  mesh_file)
virtual

This method implements reading a mesh from a specified file.

Implements libMesh::MeshInput< MeshBase >.

Definition at line 2022 of file gmv_io.C.

References _next_elem_id, _read_materials(), _read_nodes(), _read_one_cell(), _read_var(), libMesh::MeshBase::allow_renumbering(), libMesh::MeshBase::clear(), libMesh::MeshInput< MeshBase >::elems_of_dimension, libMesh::err, ierr, libMesh::MeshInput< MeshBase >::mesh(), libMesh::MeshInput< MT >::mesh(), libMesh::MeshBase::mesh_dimension(), libMesh::Trees::NODES, libMesh::MeshBase::prepare_for_use(), libMesh::processor_id(), and libMesh::MeshBase::set_mesh_dimension().

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

2023 {
2024  // This is a serial-only process for now;
2025  // the Mesh should be read on processor 0 and
2026  // broadcast later
2027  libmesh_assert_equal_to (MeshOutput<MeshBase>::mesh().processor_id(), 0);
2028 
2029  _next_elem_id = 0;
2030 
2031  libmesh_experimental();
2032 
2033 #ifndef LIBMESH_HAVE_GMV
2034 
2035  libmesh_error_msg("Cannot read GMV file " << name << " without the GMV API.");
2036 
2037 #else
2038  // We use the file-scope global variable eletypes for mapping nodes
2039  // from GMV to libmesh indices, so initialize that data now.
2040  init_eletypes();
2041 
2042  // Clear the mesh so we are sure to start from a pristeen state.
2044  mesh.clear();
2045 
2046  // Keep track of what kinds of elements this file contains
2047  elems_of_dimension.clear();
2048  elems_of_dimension.resize(4, false);
2049 
2050  // It is apparently possible for gmv files to contain
2051  // a "fromfile" directive (?) But we currently don't make
2052  // any use of this feature in LibMesh. Nonzero return val
2053  // from any function usually means an error has occurred.
2054  int ierr = GMVLib::gmvread_open_fromfileskip(const_cast<char *>(name.c_str()));
2055  if (ierr != 0)
2056  libmesh_error_msg("GMVLib::gmvread_open_fromfileskip failed!");
2057 
2058 
2059  // Loop through file until GMVEND.
2060  int iend = 0;
2061  while (iend == 0)
2062  {
2063  GMVLib::gmvread_data();
2064 
2065  // Check for GMVEND.
2066  if (GMVLib::gmv_data.keyword == GMVEND)
2067  {
2068  iend = 1;
2069  GMVLib::gmvread_close();
2070  break;
2071  }
2072 
2073  // Check for GMVERROR.
2074  if (GMVLib::gmv_data.keyword == GMVERROR)
2075  libmesh_error_msg("Encountered GMVERROR while reading!");
2076 
2077  // Process the data.
2078  switch (GMVLib::gmv_data.keyword)
2079  {
2080  case NODES:
2081  {
2082  //libMesh::out << "Reading nodes." << std::endl;
2083 
2084  if (GMVLib::gmv_data.num2 == NODES)
2085  this->_read_nodes();
2086 
2087  else if (GMVLib::gmv_data.num2 == NODE_V)
2088  libmesh_error_msg("Unsupported GMV data type NODE_V!");
2089 
2090  break;
2091  }
2092 
2093  case CELLS:
2094  {
2095  // Read 1 cell at a time
2096  // libMesh::out << "\nReading one cell." << std::endl;
2097  this->_read_one_cell();
2098  break;
2099  }
2100 
2101  case MATERIAL:
2102  {
2103  // keyword == 6
2104  // These are the materials, which we use to specify the mesh
2105  // partitioning.
2106  this->_read_materials();
2107  break;
2108  }
2109 
2110  case VARIABLE:
2111  {
2112  // keyword == 8
2113  // This is a field variable.
2114 
2115  // Check to see if we're done reading variables and break out.
2116  if (GMVLib::gmv_data.datatype == ENDKEYWORD)
2117  {
2118  // libMesh::out << "Done reading GMV variables." << std::endl;
2119  break;
2120  }
2121 
2122  if (GMVLib::gmv_data.datatype == NODE)
2123  {
2124  // libMesh::out << "Reading node field data for variable "
2125  // << GMVLib::gmv_data.name1 << std::endl;
2126  this->_read_var();
2127  break;
2128  }
2129 
2130  else
2131  {
2132  libMesh::err << "Warning: Skipping variable: "
2133  << GMVLib::gmv_data.name1
2134  << " which is of unsupported GMV datatype "
2135  << GMVLib::gmv_data.datatype
2136  << ". Nodal field data is currently the only type currently supported."
2137  << std::endl;
2138  break;
2139  }
2140 
2141  }
2142 
2143  default:
2144  libmesh_error_msg("Encountered unknown GMV keyword " << GMVLib::gmv_data.keyword);
2145 
2146  } // end switch
2147  } // end while
2148 
2149  // Set the mesh dimension to the largest encountered for an element
2150  for (unsigned char i=0; i!=4; ++i)
2151  if (elems_of_dimension[i])
2152  mesh.set_mesh_dimension(i);
2153 
2154 #if LIBMESH_DIM < 3
2155  if (mesh.mesh_dimension() > LIBMESH_DIM)
2156  libmesh_error_msg("Cannot open dimension " \
2157  << mesh.mesh_dimension() \
2158  << " mesh file when configured without " \
2159  << mesh.mesh_dimension() \
2160  << "D support.");
2161 #endif
2162 
2163  // Done reading in the mesh, now call find_neighbors, etc.
2164  // mesh.find_neighbors();
2165 
2166  // Don't allow renumbering of nodes and elements when calling
2167  // prepare_for_use(). It is usually confusing for users when the
2168  // Mesh gets renumbered automatically, since sometimes there are
2169  // other parts of the simulation which rely on the original
2170  // ordering.
2171  mesh.allow_renumbering(false);
2172  mesh.prepare_for_use();
2173 #endif
2174 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
void _read_var()
Definition: gmv_io.C:2179
void _read_nodes()
Definition: gmv_io.C:2234
void allow_renumbering(bool allow)
Definition: mesh_base.h:727
std::vector< bool > elems_of_dimension
Definition: mesh_input.h:97
unsigned int _next_elem_id
Definition: gmv_io.h:231
Base class for Mesh.
Definition: mesh_base.h:67
void prepare_for_use(const bool skip_renumber_nodes_and_elements=false, const bool skip_find_neighbors=false)
Definition: mesh_base.C:173
OStreamProxy err(std::cerr)
void set_mesh_dimension(unsigned char d)
Definition: mesh_base.h:189
void _read_materials()
Definition: gmv_io.C:2191
virtual void clear()
Definition: mesh_base.C:282
PetscErrorCode ierr
unsigned int mesh_dimension() const
Definition: mesh_base.C:147
void _read_one_cell()
Definition: gmv_io.C:2265
processor_id_type processor_id()
Definition: libmesh_base.h:96
void libMesh::MeshInput< MeshBase >::set_n_partitions ( unsigned int  n_parts)
inlineprotectedinherited

Sets the number of partitions in the mesh. Typically this gets done by the partitioner, but some parallel file formats begin "pre-partitioned".

Definition at line 91 of file mesh_input.h.

References libMesh::MeshInput< MT >::mesh().

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

91 { this->mesh().set_n_partitions() = n_parts; }
unsigned int & set_n_partitions()
Definition: mesh_base.h:1295
void libMesh::MeshInput< MeshBase >::skip_comment_lines ( std::istream &  in,
const char  comment_start 
)
protectedinherited

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

bool& libMesh::GMVIO::subdivide_second_order ( )
inline

Flag indicating whether or not to subdivide second order elements

Definition at line 120 of file gmv_io.h.

References _subdivide_second_order.

Referenced by write_ascii_new_impl(), and write_ascii_old_impl().

120 { return _subdivide_second_order; }
bool _subdivide_second_order
Definition: gmv_io.h:212
void libMesh::GMVIO::write ( const std::string &  fname)
virtual

This method implements writing a mesh to a specified file.

Implements libMesh::MeshOutput< MeshBase >.

Definition at line 267 of file gmv_io.C.

References binary(), write_ascii_old_impl(), and write_binary().

Referenced by libMesh::NameBasedIO::write().

268 {
269  if (this->binary())
270  this->write_binary (fname);
271  else
272  this->write_ascii_old_impl (fname);
273 }
bool & binary()
Definition: gmv_io.h:95
void write_ascii_old_impl(const std::string &, const std::vector< Number > *=libmesh_nullptr, const std::vector< std::string > *=libmesh_nullptr)
Definition: gmv_io.C:587
void write_binary(const std::string &, const std::vector< Number > *=libmesh_nullptr, const std::vector< std::string > *=libmesh_nullptr)
Definition: gmv_io.C:1275
void libMesh::GMVIO::write_ascii_new_impl ( const std::string &  fname,
const std::vector< Number > *  v = libmesh_nullptr,
const std::vector< std::string > *  solution_names = libmesh_nullptr 
)

This method implements writing a mesh with nodal data to a specified file where the nodal data and variable names are optionally provided. This will write an ASCII file. This is the new implementation (without subcells).

Definition at line 291 of file gmv_io.C.

References _cell_centered_data, std::abs(), libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::MeshOutput< MeshBase >::ascii_precision(), end, libMesh::err, libMesh::DofObject::id(), libMesh::libmesh_assert(), libmesh_nullptr, std::max(), libMesh::MeshInput< MeshBase >::mesh(), libMesh::MeshOutput< MT >::mesh(), libMesh::MeshBase::n_active_elem(), libMesh::MeshBase::n_nodes(), libMesh::Elem::n_nodes(), libMesh::MeshBase::n_partitions(), libMesh::Elem::n_sub_elem(), n_vars, libMesh::Elem::node_id(), libMesh::Elem::p_level(), p_levels(), partitioning(), libMesh::MeshBase::point(), libMesh::processor_id(), libMesh::Real, subdivide_second_order(), libMesh::Elem::type(), write_ascii_old_impl(), and write_subdomain_id_as_material().

Referenced by p_levels().

294 {
295 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
296 
297  libMesh::err << "WARNING: GMVIO::write_ascii_new_impl() not infinite-element aware!"
298  << std::endl;
299  libmesh_here();
300 
301  // Set it to our current precision
302  this->write_ascii_old_impl (fname, v, solution_names);
303 
304 #else
305 
306  // Get a reference to the mesh
308 
309  // This is a parallel_only function
310  const unsigned int n_active_elem = mesh.n_active_elem();
311 
313  return;
314 
315  // Open the output file stream
316  std::ofstream out_stream (fname.c_str());
317 
318  out_stream << std::setprecision(this->ascii_precision());
319 
320  // Make sure it opened correctly
321  if (!out_stream.good())
322  libmesh_file_error(fname.c_str());
323 
324  unsigned int mesh_max_p_level = 0;
325 
326  // Begin interfacing with the GMV data file
327  {
328  out_stream << "gmvinput ascii\n\n";
329 
330  // write the nodes
331  out_stream << "nodes " << mesh.n_nodes() << "\n";
332  for (unsigned int n=0; n<mesh.n_nodes(); n++)
333  out_stream << mesh.point(n)(0) << " ";
334  out_stream << "\n";
335 
336  for (unsigned int n=0; n<mesh.n_nodes(); n++)
337 #if LIBMESH_DIM > 1
338  out_stream << mesh.point(n)(1) << " ";
339 #else
340  out_stream << 0. << " ";
341 #endif
342  out_stream << "\n";
343 
344  for (unsigned int n=0; n<mesh.n_nodes(); n++)
345 #if LIBMESH_DIM > 2
346  out_stream << mesh.point(n)(2) << " ";
347 #else
348  out_stream << 0. << " ";
349 #endif
350  out_stream << "\n\n";
351  }
352 
353  {
354  // write the connectivity
355  out_stream << "cells " << n_active_elem << "\n";
356 
357  // initialize the eletypes map (eletypes is a file-global variable)
358  init_eletypes();
359 
362 
363  for ( ; it != end; ++it)
364  {
365  const Elem * elem = *it;
366 
367  mesh_max_p_level = std::max(mesh_max_p_level,
368  elem->p_level());
369 
370  // Make sure we have a valid entry for
371  // the current element type.
372  libmesh_assert (eletypes.count(elem->type()));
373 
374  const ElementDefinition & ele = eletypes[elem->type()];
375 
376  // The element mapper better not require any more nodes
377  // than are present in the current element!
378  libmesh_assert_less_equal (ele.node_map.size(), elem->n_nodes());
379 
380  out_stream << ele.label << "\n";
381  for (std::size_t i=0; i < ele.node_map.size(); i++)
382  out_stream << elem->node_id(ele.node_map[i])+1 << " ";
383  out_stream << "\n";
384  }
385  out_stream << "\n";
386  }
387 
388  // optionally write the partition information
389  if (this->partitioning())
390  {
391  if (this->write_subdomain_id_as_material())
392  libmesh_error_msg("Not yet supported in GMVIO::write_ascii_new_impl");
393 
394  else // write processor IDs as materials. This is the default
395  {
396  out_stream << "material "
397  << mesh.n_partitions()
398  // Note: GMV may give you errors like
399  // Error, material for cell 1 is greater than 1
400  // Error, material for cell 2 is greater than 1
401  // Error, material for cell 3 is greater than 1
402  // ... because you put the wrong number of partitions here.
403  // To ensure you write the correct number of materials, call
404  // mesh.recalculate_n_partitions() before writing out the
405  // mesh.
406  // Note: we can't call it now because the Mesh is const here and
407  // it is a non-const function.
408  << " 0\n";
409 
410  for (unsigned int proc=0; proc<mesh.n_partitions(); proc++)
411  out_stream << "proc_" << proc << "\n";
412 
415 
416  // FIXME - don't we need to use an ElementDefinition here? - RHS
417  for ( ; it != end; ++it)
418  out_stream << (*it)->processor_id()+1 << "\n";
419  out_stream << "\n";
420  }
421  }
422 
423  // If there are *any* variables at all in the system (including
424  // p level, or arbitrary cell-based data)
425  // to write, the gmv file needs to contain the word "variable"
426  // on a line by itself.
427  bool write_variable = false;
428 
429  // 1.) p-levels
430  if (this->p_levels() && mesh_max_p_level)
431  write_variable = true;
432 
433  // 2.) solution data
434  if ((solution_names != libmesh_nullptr) && (v != libmesh_nullptr))
435  write_variable = true;
436 
437  // 3.) cell-centered data
438  if (!(this->_cell_centered_data.empty()))
439  write_variable = true;
440 
441  if (write_variable)
442  out_stream << "variable\n";
443 
444  // if ((this->p_levels() && mesh_max_p_level) ||
445  // ((solution_names != libmesh_nullptr) && (v != libmesh_nullptr)))
446  // out_stream << "variable\n";
447 
448  // optionally write the polynomial degree information
449  if (this->p_levels() && mesh_max_p_level)
450  {
451  out_stream << "p_level 0\n";
452 
455 
456  for ( ; it != end; ++it)
457  {
458  const Elem * elem = *it;
459 
460  const ElementDefinition & ele = eletypes[elem->type()];
461 
462  // The element mapper better not require any more nodes
463  // than are present in the current element!
464  libmesh_assert_less_equal (ele.node_map.size(), elem->n_nodes());
465 
466  for (std::size_t i=0; i < ele.node_map.size(); i++)
467  out_stream << elem->p_level() << " ";
468  }
469  out_stream << "\n\n";
470  }
471 
472 
473  // optionally write cell-centered data
474  if (!(this->_cell_centered_data.empty()))
475  {
476  std::map<std::string, const std::vector<Real> *>::iterator it = this->_cell_centered_data.begin();
477  const std::map<std::string, const std::vector<Real> *>::iterator end = this->_cell_centered_data.end();
478 
479  for (; it != end; ++it)
480  {
481  // write out the variable name, followed by a zero.
482  out_stream << (*it).first << " 0\n";
483 
484  const std::vector<Real> * the_array = (*it).second;
485 
486  // Loop over active elements, write out cell data. If second-order cells
487  // are split into sub-elements, the sub-elements inherit their parent's
488  // cell-centered data.
491 
492  for (; elem_it != elem_end; ++elem_it)
493  {
494  const Elem * e = *elem_it;
495 
496  // Use the element's ID to find the value.
497  libmesh_assert_less (e->id(), the_array->size());
498  const Real the_value = the_array->operator[](e->id());
499 
500  if (this->subdivide_second_order())
501  for (unsigned int se=0; se < e->n_sub_elem(); se++)
502  out_stream << the_value << " ";
503  else
504  out_stream << the_value << " ";
505  }
506 
507  out_stream << "\n\n";
508  }
509  }
510 
511 
512  // optionally write the data
513  if ((solution_names != libmesh_nullptr) && (v != libmesh_nullptr))
514  {
515  const unsigned int n_vars = solution_names->size();
516 
517  if (!(v->size() == mesh.n_nodes()*n_vars))
518  libMesh::err << "ERROR: v->size()=" << v->size()
519  << ", mesh.n_nodes()=" << mesh.n_nodes()
520  << ", n_vars=" << n_vars
521  << ", mesh.n_nodes()*n_vars=" << mesh.n_nodes()*n_vars
522  << "\n";
523 
524  libmesh_assert_equal_to (v->size(), mesh.n_nodes()*n_vars);
525 
526  for (unsigned int c=0; c<n_vars; c++)
527  {
528 
529 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
530 
531  // in case of complex data, write _three_ data sets
532  // for each component
533 
534  // this is the real part
535  out_stream << "r_" << (*solution_names)[c] << " 1\n";
536 
537  for (unsigned int n=0; n<mesh.n_nodes(); n++)
538  out_stream << (*v)[n*n_vars + c].real() << " ";
539 
540  out_stream << "\n\n";
541 
542  // this is the imaginary part
543  out_stream << "i_" << (*solution_names)[c] << " 1\n";
544 
545  for (unsigned int n=0; n<mesh.n_nodes(); n++)
546  out_stream << (*v)[n*n_vars + c].imag() << " ";
547 
548  out_stream << "\n\n";
549 
550  // this is the magnitude
551  out_stream << "a_" << (*solution_names)[c] << " 1\n";
552  for (unsigned int n=0; n<mesh.n_nodes(); n++)
553  out_stream << std::abs((*v)[n*n_vars + c]) << " ";
554 
555  out_stream << "\n\n";
556 
557 #else
558 
559  out_stream << (*solution_names)[c] << " 1\n";
560 
561  for (unsigned int n=0; n<mesh.n_nodes(); n++)
562  out_stream << (*v)[n*n_vars + c] << " ";
563 
564  out_stream << "\n\n";
565 
566 #endif
567  }
568 
569  }
570 
571  // If we wrote any variables, we have to close the variable section now
572  if (write_variable)
573  out_stream << "endvars\n";
574 
575 
576  // end of the file
577  out_stream << "\nendgmv\n";
578 
579 #endif
580 }
double abs(double a)
virtual dof_id_type n_active_elem() const =0
virtual const Point & point(const dof_id_type i) const =0
virtual ElemType type() const =0
unsigned int p_level() const
Definition: elem.h:2218
unsigned int n_partitions() const
Definition: mesh_base.h:811
The base class for all geometric element types.
Definition: elem.h:86
const class libmesh_nullptr_t libmesh_nullptr
unsigned int & ascii_precision()
IterBase * end
const MT & mesh() const
Definition: mesh_output.h:216
const unsigned int n_vars
Definition: tecplot_io.C:68
bool & partitioning()
Definition: gmv_io.h:107
long double max(long double a, double b)
Base class for Mesh.
Definition: mesh_base.h:67
libmesh_assert(j)
virtual unsigned int n_nodes() const =0
void write_ascii_old_impl(const std::string &, const std::vector< Number > *=libmesh_nullptr, const std::vector< std::string > *=libmesh_nullptr)
Definition: gmv_io.C:587
bool & write_subdomain_id_as_material()
Definition: gmv_io.h:114
std::map< std::string, const std::vector< Real > * > _cell_centered_data
Definition: gmv_io.h:225
bool & p_levels()
Definition: gmv_io.h:126
OStreamProxy err(std::cerr)
virtual element_iterator active_elements_begin()=0
virtual element_iterator active_elements_end()=0
bool & subdivide_second_order()
Definition: gmv_io.h:120
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
dof_id_type node_id(const unsigned int i) const
Definition: elem.h:1683
virtual unsigned int n_sub_elem() const =0
dof_id_type id() const
Definition: dof_object.h:624
virtual dof_id_type n_nodes() const =0
processor_id_type processor_id()
Definition: libmesh_base.h:96
void libMesh::GMVIO::write_ascii_old_impl ( const std::string &  fname,
const std::vector< Number > *  v = libmesh_nullptr,
const std::vector< std::string > *  solution_names = libmesh_nullptr 
)
private

This method implements writing a mesh with nodal data to a specified file where the nodal data and variable names are optionally provided. This will write an ASCII file. This is the old implementation (using subcells) which was the default in libMesh-0.4.3-rc2.

Definition at line 587 of file gmv_io.C.

References _cell_centered_data, std::abs(), libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::MeshOutput< MeshBase >::ascii_precision(), libMesh::Elem::build(), libMesh::Elem::dim(), end, libMesh::Utility::enum_to_string(), libMesh::err, libMesh::FIRST, libMesh::Elem::first_order_equivalent_type(), libMesh::HEX20, libMesh::HEX27, libMesh::HEX8, libMesh::DofObject::id(), libMesh::INFHEX16, libMesh::INFHEX18, libMesh::INFHEX8, libMesh::INFPRISM12, libMesh::INFPRISM6, libMesh::INFQUAD4, libMesh::INFQUAD6, libmesh_nullptr, std::max(), libMesh::MeshBase::max_node_id(), libMesh::MeshInput< MeshBase >::mesh(), libMesh::MeshOutput< MT >::mesh(), libMesh::MeshBase::n_active_elem(), libMesh::MeshBase::n_active_sub_elem(), libMesh::MeshBase::n_nodes(), libMesh::MeshBase::n_partitions(), libMesh::Elem::n_sub_elem(), n_vars, libMesh::Elem::p_level(), p_levels(), partitioning(), libMesh::MeshBase::point(), libMesh::PRISM15, libMesh::PRISM18, libMesh::PRISM6, libMesh::processor_id(), libMesh::PYRAMID5, libMesh::QUAD4, libMesh::QUAD8, libMesh::QUAD9, libMesh::Real, subdivide_second_order(), libMesh::TECPLOT, libMesh::TET10, libMesh::TET4, libMesh::TRI3, libMesh::TRI6, and write_subdomain_id_as_material().

Referenced by p_levels(), write(), write_ascii_new_impl(), and write_nodal_data().

590 {
591  // Get a reference to the mesh
593 
594  // Use a MeshSerializer object to gather a parallel mesh before outputting it.
595  // Note that we cast away constness here (which is bad), but the destructor of
596  // the MeshSerializer object reparallelizes the Mesh, hopefully keeping it
597  // "logically const" outside the context of this function...
598  MeshSerializer serialize(const_cast<MeshBase &>(mesh),
600 
601  // These are parallel_only functions
602  const dof_id_type n_active_elem = mesh.n_active_elem(),
603  n_active_sub_elem = mesh.n_active_sub_elem();
604 
606  return;
607 
608  // Open the output file stream
609  std::ofstream out_stream (fname.c_str());
610 
611  // Set it to our current precision
612  out_stream << std::setprecision(this->ascii_precision());
613 
614  // Make sure it opened correctly
615  if (!out_stream.good())
616  libmesh_file_error(fname.c_str());
617 
618  // Make sure our nodes are contiguous and serialized
619  libmesh_assert_equal_to (mesh.n_nodes(), mesh.max_node_id());
620 
621  // libmesh_assert (mesh.is_serial());
622  // if (!mesh.is_serial())
623  // {
624  // if (MeshOutput<MeshBase>::mesh().processor_id() == 0)
625  // libMesh::err << "Error: GMVIO cannot yet write a DistributedMesh solution"
626  // << std::endl;
627  // return;
628  // }
629 
630  unsigned int mesh_max_p_level = 0;
631 
632  // Begin interfacing with the GMV data file
633 
634  // FIXME - if subdivide_second_order() is off,
635  // we probably should only be writing the
636  // vertex nodes - RHS
637  {
638  // write the nodes
639 
640  out_stream << "gmvinput ascii\n\n";
641  out_stream << "nodes " << mesh.n_nodes() << '\n';
642  for (unsigned int n=0; n<mesh.n_nodes(); n++)
643  out_stream << mesh.point(n)(0) << " ";
644 
645  out_stream << '\n';
646 
647  for (unsigned int n=0; n<mesh.n_nodes(); n++)
648 #if LIBMESH_DIM > 1
649  out_stream << mesh.point(n)(1) << " ";
650 #else
651  out_stream << 0. << " ";
652 #endif
653 
654  out_stream << '\n';
655 
656  for (unsigned int n=0; n<mesh.n_nodes(); n++)
657 #if LIBMESH_DIM > 2
658  out_stream << mesh.point(n)(2) << " ";
659 #else
660  out_stream << 0. << " ";
661 #endif
662 
663  out_stream << '\n' << '\n';
664  }
665 
666 
667 
668  {
669  // write the connectivity
670 
671  out_stream << "cells ";
672  if (this->subdivide_second_order())
673  out_stream << n_active_sub_elem;
674  else
675  out_stream << n_active_elem;
676  out_stream << '\n';
677 
678  // The same temporary storage will be used for each element
679  std::vector<dof_id_type> conn;
680 
683 
684  for ( ; it != end; ++it)
685  {
686  const Elem * elem = *it;
687 
688  mesh_max_p_level = std::max(mesh_max_p_level,
689  elem->p_level());
690 
691  switch (elem->dim())
692  {
693  case 1:
694  {
695  if (this->subdivide_second_order())
696  for (unsigned int se=0; se<(*it)->n_sub_elem(); se++)
697  {
698  out_stream << "line 2\n";
699  (*it)->connectivity(se, TECPLOT, conn);
700  for (std::size_t i=0; i<conn.size(); i++)
701  out_stream << conn[i] << " ";
702 
703  out_stream << '\n';
704  }
705  else
706  {
707  out_stream << "line 2\n";
708  if ((*it)->default_order() == FIRST)
709  (*it)->connectivity(0, TECPLOT, conn);
710  else
711  {
713  for (unsigned int i = 0; i != lo_elem->n_nodes(); ++i)
714  lo_elem->set_node(i) = (*it)->node_ptr(i);
715  lo_elem->connectivity(0, TECPLOT, conn);
716  }
717  for (std::size_t i=0; i<conn.size(); i++)
718  out_stream << conn[i] << " ";
719 
720  out_stream << '\n';
721  }
722 
723  break;
724  }
725 
726  case 2:
727  {
728  if (this->subdivide_second_order())
729  for (unsigned int se=0; se<(*it)->n_sub_elem(); se++)
730  {
731  // Quad elements
732  if (((*it)->type() == QUAD4) ||
733  ((*it)->type() == QUAD8) || // Note: QUAD8 will be output as one central quad and
734  // four surrounding triangles (though they will be written
735  // to GMV as QUAD4s).
736  ((*it)->type() == QUAD9)
737 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
738  || ((*it)->type() == INFQUAD4)
739  || ((*it)->type() == INFQUAD6)
740 #endif
741  )
742  {
743  out_stream << "quad 4\n";
744  (*it)->connectivity(se, TECPLOT, conn);
745  for (std::size_t i=0; i<conn.size(); i++)
746  out_stream << conn[i] << " ";
747  }
748 
749  // Triangle elements
750  else if (((*it)->type() == TRI3) ||
751  ((*it)->type() == TRI6))
752  {
753  out_stream << "tri 3\n";
754  (*it)->connectivity(se, TECPLOT, conn);
755  for (unsigned int i=0; i<3; i++)
756  out_stream << conn[i] << " ";
757  }
758  else
759  libmesh_error_msg("Unsupported element type: " << Utility::enum_to_string((*it)->type()));
760  }
761  else // !this->subdivide_second_order()
762  {
763  // Quad elements
764  if (((*it)->type() == QUAD4)
765 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
766  || ((*it)->type() == INFQUAD4)
767 #endif
768  )
769  {
770  (*it)->connectivity(0, TECPLOT, conn);
771  out_stream << "quad 4\n";
772  for (std::size_t i=0; i<conn.size(); i++)
773  out_stream << conn[i] << " ";
774  }
775  else if (((*it)->type() == QUAD8) ||
776  ((*it)->type() == QUAD9)
777 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
778  || ((*it)->type() == INFQUAD6)
779 #endif
780  )
781  {
783  for (unsigned int i = 0; i != lo_elem->n_nodes(); ++i)
784  lo_elem->set_node(i) = (*it)->node_ptr(i);
785  lo_elem->connectivity(0, TECPLOT, conn);
786  out_stream << "quad 4\n";
787  for (std::size_t i=0; i<conn.size(); i++)
788  out_stream << conn[i] << " ";
789  }
790  else if ((*it)->type() == TRI3)
791  {
792  (*it)->connectivity(0, TECPLOT, conn);
793  out_stream << "tri 3\n";
794  for (unsigned int i=0; i<3; i++)
795  out_stream << conn[i] << " ";
796  }
797  else if ((*it)->type() == TRI6)
798  {
800  for (unsigned int i = 0; i != lo_elem->n_nodes(); ++i)
801  lo_elem->set_node(i) = (*it)->node_ptr(i);
802  lo_elem->connectivity(0, TECPLOT, conn);
803  out_stream << "tri 3\n";
804  for (unsigned int i=0; i<3; i++)
805  out_stream << conn[i] << " ";
806  }
807 
808  out_stream << '\n';
809  }
810 
811  break;
812  }
813 
814  case 3:
815  {
816  if (this->subdivide_second_order())
817  for (unsigned int se=0; se<(*it)->n_sub_elem(); se++)
818  {
819 
820 #ifndef LIBMESH_ENABLE_INFINITE_ELEMENTS
821  if (((*it)->type() == HEX8) ||
822  ((*it)->type() == HEX27))
823  {
824  out_stream << "phex8 8\n";
825  (*it)->connectivity(se, TECPLOT, conn);
826  for (std::size_t i=0; i<conn.size(); i++)
827  out_stream << conn[i] << " ";
828  }
829 
830  else if ((*it)->type() == HEX20)
831  {
832  out_stream << "phex20 20\n";
833  out_stream << (*it)->node_id(0)+1 << " "
834  << (*it)->node_id(1)+1 << " "
835  << (*it)->node_id(2)+1 << " "
836  << (*it)->node_id(3)+1 << " "
837  << (*it)->node_id(4)+1 << " "
838  << (*it)->node_id(5)+1 << " "
839  << (*it)->node_id(6)+1 << " "
840  << (*it)->node_id(7)+1 << " "
841  << (*it)->node_id(8)+1 << " "
842  << (*it)->node_id(9)+1 << " "
843  << (*it)->node_id(10)+1 << " "
844  << (*it)->node_id(11)+1 << " "
845  << (*it)->node_id(16)+1 << " "
846  << (*it)->node_id(17)+1 << " "
847  << (*it)->node_id(18)+1 << " "
848  << (*it)->node_id(19)+1 << " "
849  << (*it)->node_id(12)+1 << " "
850  << (*it)->node_id(13)+1 << " "
851  << (*it)->node_id(14)+1 << " "
852  << (*it)->node_id(15)+1 << " ";
853  }
854 
855  // According to our copy of gmvread.c, this is the
856  // mapping for the Hex27 element. Unfortunately,
857  // I tried it and Paraview does not seem to be
858  // able to read Hex27 elements. Since this is
859  // unlikely to change any time soon, we'll
860  // continue to write out Hex27 elements as 8 Hex8
861  // sub-elements.
862  //
863  // TODO:
864  // 1.) If we really wanted to use this code for
865  // something, we'd want to avoid repeating the
866  // hard-coded node ordering from the HEX20 case.
867  // These should both be able to use
868  // ElementDefinitions.
869  // 2.) You would also need to change
870  // Hex27::n_sub_elem() to return 1, not sure how
871  // much other code that would affect...
872 
873  // else if ((*it)->type() == HEX27)
874  // {
875  // out_stream << "phex27 27\n";
876  // out_stream << (*it)->node_id(0)+1 << " "
877  // << (*it)->node_id(1)+1 << " "
878  // << (*it)->node_id(2)+1 << " "
879  // << (*it)->node_id(3)+1 << " "
880  // << (*it)->node_id(4)+1 << " "
881  // << (*it)->node_id(5)+1 << " "
882  // << (*it)->node_id(6)+1 << " "
883  // << (*it)->node_id(7)+1 << " "
884  // << (*it)->node_id(8)+1 << " "
885  // << (*it)->node_id(9)+1 << " "
886  // << (*it)->node_id(10)+1 << " "
887  // << (*it)->node_id(11)+1 << " "
888  // << (*it)->node_id(16)+1 << " "
889  // << (*it)->node_id(17)+1 << " "
890  // << (*it)->node_id(18)+1 << " "
891  // << (*it)->node_id(19)+1 << " "
892  // << (*it)->node_id(12)+1 << " "
893  // << (*it)->node_id(13)+1 << " "
894  // << (*it)->node_id(14)+1 << " "
895  // << (*it)->node_id(15)+1 << " "
896  // // mid-face nodes
897  // << (*it)->node_id(21)+1 << " " // GMV front
898  // << (*it)->node_id(22)+1 << " " // GMV right
899  // << (*it)->node_id(23)+1 << " " // GMV back
900  // << (*it)->node_id(24)+1 << " " // GMV left
901  // << (*it)->node_id(20)+1 << " " // GMV bottom
902  // << (*it)->node_id(25)+1 << " " // GMV top
903  // // center node
904  // << (*it)->node_id(26)+1 << " ";
905  // }
906 
907 #else // LIBMESH_ENABLE_INFINITE_ELEMENTS
908 
909  // In case of infinite elements, HEX20
910  // should be handled just like the
911  // INFHEX16, since these connect to each other
912  if (((*it)->type() == HEX8) ||
913  ((*it)->type() == HEX27) ||
914  ((*it)->type() == INFHEX8) ||
915  ((*it)->type() == INFHEX16) ||
916  ((*it)->type() == INFHEX18) ||
917  ((*it)->type() == HEX20))
918  {
919  out_stream << "phex8 8\n";
920  (*it)->connectivity(se, TECPLOT, conn);
921  for (std::size_t i=0; i<conn.size(); i++)
922  out_stream << conn[i] << " ";
923  }
924 #endif
925 
926  else if (((*it)->type() == TET4) ||
927  ((*it)->type() == TET10))
928  {
929  out_stream << "tet 4\n";
930  // Tecplot connectivity returns 8 entries for
931  // the Tet, enough to store it as a degenerate Hex.
932  // For GMV we only pick out the four relevant node
933  // indices.
934  (*it)->connectivity(se, TECPLOT, conn);
935  out_stream << conn[0] << " " // libmesh tet node 0
936  << conn[2] << " " // libmesh tet node 2
937  << conn[1] << " " // libmesh tet node 1
938  << conn[4] << " "; // libmesh tet node 3
939  }
940 #ifndef LIBMESH_ENABLE_INFINITE_ELEMENTS
941  else if (((*it)->type() == PRISM6) ||
942  ((*it)->type() == PRISM15) ||
943  ((*it)->type() == PRISM18) ||
944  ((*it)->type() == PYRAMID5))
945 #else
946  else if (((*it)->type() == PRISM6) ||
947  ((*it)->type() == PRISM15) ||
948  ((*it)->type() == PRISM18) ||
949  ((*it)->type() == PYRAMID5) ||
950  ((*it)->type() == INFPRISM6) ||
951  ((*it)->type() == INFPRISM12))
952 #endif
953  {
954  // Note that the prisms are treated as
955  // degenerated phex8's.
956  out_stream << "phex8 8\n";
957  (*it)->connectivity(se, TECPLOT, conn);
958  for (std::size_t i=0; i<conn.size(); i++)
959  out_stream << conn[i] << " ";
960  }
961 
962  else
963  libmesh_error_msg("Encountered an unrecognized element " \
964  << "type: " << (*it)->type() \
965  << "\nPossibly a dim-1 dimensional " \
966  << "element? Aborting...");
967 
968  out_stream << '\n';
969  }
970  else // !this->subdivide_second_order()
971  {
973  for (unsigned int i = 0; i != lo_elem->n_nodes(); ++i)
974  lo_elem->set_node(i) = (*it)->node_ptr(i);
975  if ((lo_elem->type() == HEX8)
976 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
977  || (lo_elem->type() == HEX27)
978 #endif
979  )
980  {
981  out_stream << "phex8 8\n";
982  lo_elem->connectivity(0, TECPLOT, conn);
983  for (std::size_t i=0; i<conn.size(); i++)
984  out_stream << conn[i] << " ";
985  }
986 
987  else if (lo_elem->type() == TET4)
988  {
989  out_stream << "tet 4\n";
990  lo_elem->connectivity(0, TECPLOT, conn);
991  out_stream << conn[0] << " "
992  << conn[2] << " "
993  << conn[1] << " "
994  << conn[4] << " ";
995  }
996  else if ((lo_elem->type() == PRISM6)
997 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
998  || (lo_elem->type() == INFPRISM6)
999 #endif
1000  )
1001  {
1002  // Note that the prisms are treated as
1003  // degenerated phex8's.
1004  out_stream << "phex8 8\n";
1005  lo_elem->connectivity(0, TECPLOT, conn);
1006  for (std::size_t i=0; i<conn.size(); i++)
1007  out_stream << conn[i] << " ";
1008  }
1009 
1010  else
1011  libmesh_error_msg("Encountered an unrecognized element " \
1012  << "type. Possibly a dim-1 dimensional " \
1013  << "element? Aborting...");
1014 
1015  out_stream << '\n';
1016  }
1017 
1018  break;
1019  }
1020 
1021  default:
1022  libmesh_error_msg("Unsupported element dimension: " <<
1023  elem->dim());
1024  }
1025  }
1026 
1027  out_stream << '\n';
1028  }
1029 
1030 
1031 
1032  // optionally write the partition information
1033  if (this->partitioning())
1034  {
1035  if (this->write_subdomain_id_as_material())
1036  {
1037  // Subdomain IDs can be non-contiguous and need not
1038  // necessarily start at 0. Furthermore, since the user is
1039  // free to define subdomain_id_type to be a signed type, we
1040  // can't even assume max(subdomain_id) >= # unique subdomain ids.
1041 
1042  // We build a map<subdomain_id, unsigned> to associate to each
1043  // user-selected subdomain ID a unique, contiguous unsigned value
1044  // which we can write to file.
1045  std::map<subdomain_id_type, unsigned> sbdid_map;
1046  typedef std::map<subdomain_id_type, unsigned>::iterator sbdid_map_iter;
1047  {
1050 
1051  for ( ; it != end; ++it)
1052  {
1053  // Try to insert with dummy value
1054  sbdid_map.insert( std::make_pair((*it)->subdomain_id(), 0) );
1055  }
1056  }
1057 
1058  // Map is created, iterate through it to set indices. They will be
1059  // used repeatedly below.
1060  {
1061  unsigned ctr=0;
1062  for (sbdid_map_iter it=sbdid_map.begin(); it != sbdid_map.end(); ++it)
1063  (*it).second = ctr++;
1064  }
1065 
1066  out_stream << "material "
1067  << sbdid_map.size()
1068  << " 0\n";
1069 
1070  for (std::size_t sbdid=0; sbdid<sbdid_map.size(); sbdid++)
1071  out_stream << "proc_" << sbdid << "\n";
1072 
1075 
1076  for ( ; it != end; ++it)
1077  {
1078  // Find the unique index for (*it)->subdomain_id(), print that to file
1079  sbdid_map_iter map_iter = sbdid_map.find( (*it)->subdomain_id() );
1080  unsigned gmv_mat_number = (*map_iter).second;
1081 
1082  if (this->subdivide_second_order())
1083  for (unsigned int se=0; se<(*it)->n_sub_elem(); se++)
1084  out_stream << gmv_mat_number+1 << '\n';
1085  else
1086  out_stream << gmv_mat_number+1 << "\n";
1087  }
1088  out_stream << '\n';
1089 
1090  }
1091  else // write processor IDs as materials. This is the default
1092  {
1093  out_stream << "material "
1094  << mesh.n_partitions()
1095  << " 0"<< '\n';
1096 
1097  for (unsigned int proc=0; proc<mesh.n_partitions(); proc++)
1098  out_stream << "proc_" << proc << '\n';
1099 
1102 
1103  for ( ; it != end; ++it)
1104  if (this->subdivide_second_order())
1105  for (unsigned int se=0; se<(*it)->n_sub_elem(); se++)
1106  out_stream << (*it)->processor_id()+1 << '\n';
1107  else
1108  out_stream << (*it)->processor_id()+1 << '\n';
1109 
1110  out_stream << '\n';
1111  }
1112  }
1113 
1114 
1115  // If there are *any* variables at all in the system (including
1116  // p level, or arbitrary cell-based data)
1117  // to write, the gmv file needs to contain the word "variable"
1118  // on a line by itself.
1119  bool write_variable = false;
1120 
1121  // 1.) p-levels
1122  if (this->p_levels() && mesh_max_p_level)
1123  write_variable = true;
1124 
1125  // 2.) solution data
1126  if ((solution_names != libmesh_nullptr) && (v != libmesh_nullptr))
1127  write_variable = true;
1128 
1129  // 3.) cell-centered data
1130  if (!(this->_cell_centered_data.empty()))
1131  write_variable = true;
1132 
1133  if (write_variable)
1134  out_stream << "variable\n";
1135 
1136 
1137  // optionally write the p-level information
1138  if (this->p_levels() && mesh_max_p_level)
1139  {
1140  out_stream << "p_level 0\n";
1141 
1144 
1145  for ( ; it != end; ++it)
1146  if (this->subdivide_second_order())
1147  for (unsigned int se=0; se<(*it)->n_sub_elem(); se++)
1148  out_stream << (*it)->p_level() << " ";
1149  else
1150  out_stream << (*it)->p_level() << " ";
1151  out_stream << "\n\n";
1152  }
1153 
1154 
1155 
1156 
1157  // optionally write cell-centered data
1158  if (!(this->_cell_centered_data.empty()))
1159  {
1160  std::map<std::string, const std::vector<Real> *>::iterator it = this->_cell_centered_data.begin();
1161  const std::map<std::string, const std::vector<Real> *>::iterator end = this->_cell_centered_data.end();
1162 
1163  for (; it != end; ++it)
1164  {
1165  // write out the variable name, followed by a zero.
1166  out_stream << (*it).first << " 0\n";
1167 
1168  const std::vector<Real> * the_array = (*it).second;
1169 
1170  // Loop over active elements, write out cell data. If second-order cells
1171  // are split into sub-elements, the sub-elements inherit their parent's
1172  // cell-centered data.
1174  const MeshBase::const_element_iterator elem_end = mesh.active_elements_end();
1175 
1176  for (; elem_it != elem_end; ++elem_it)
1177  {
1178  const Elem * e = *elem_it;
1179 
1180  // Use the element's ID to find the value...
1181  libmesh_assert_less (e->id(), the_array->size());
1182  const Real the_value = the_array->operator[](e->id());
1183 
1184  if (this->subdivide_second_order())
1185  for (unsigned int se=0; se < e->n_sub_elem(); se++)
1186  out_stream << the_value << " ";
1187  else
1188  out_stream << the_value << " ";
1189  }
1190 
1191  out_stream << "\n\n";
1192  }
1193  }
1194 
1195 
1196 
1197 
1198  // optionally write the data
1199  if ((solution_names != libmesh_nullptr) &&
1200  (v != libmesh_nullptr))
1201  {
1202  const unsigned int n_vars =
1203  cast_int<unsigned int>(solution_names->size());
1204 
1205  if (!(v->size() == mesh.n_nodes()*n_vars))
1206  libMesh::err << "ERROR: v->size()=" << v->size()
1207  << ", mesh.n_nodes()=" << mesh.n_nodes()
1208  << ", n_vars=" << n_vars
1209  << ", mesh.n_nodes()*n_vars=" << mesh.n_nodes()*n_vars
1210  << std::endl;
1211 
1212  libmesh_assert_equal_to (v->size(), mesh.n_nodes()*n_vars);
1213 
1214  for (unsigned int c=0; c<n_vars; c++)
1215  {
1216 
1217 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
1218 
1219  // in case of complex data, write _tree_ data sets
1220  // for each component
1221 
1222  // this is the real part
1223  out_stream << "r_" << (*solution_names)[c] << " 1\n";
1224 
1225  for (unsigned int n=0; n<mesh.n_nodes(); n++)
1226  out_stream << (*v)[n*n_vars + c].real() << " ";
1227 
1228  out_stream << '\n' << '\n';
1229 
1230 
1231  // this is the imaginary part
1232  out_stream << "i_" << (*solution_names)[c] << " 1\n";
1233 
1234  for (unsigned int n=0; n<mesh.n_nodes(); n++)
1235  out_stream << (*v)[n*n_vars + c].imag() << " ";
1236 
1237  out_stream << '\n' << '\n';
1238 
1239  // this is the magnitude
1240  out_stream << "a_" << (*solution_names)[c] << " 1\n";
1241  for (unsigned int n=0; n<mesh.n_nodes(); n++)
1242  out_stream << std::abs((*v)[n*n_vars + c]) << " ";
1243 
1244  out_stream << '\n' << '\n';
1245 
1246 #else
1247 
1248  out_stream << (*solution_names)[c] << " 1\n";
1249 
1250  for (unsigned int n=0; n<mesh.n_nodes(); n++)
1251  out_stream << (*v)[n*n_vars + c] << " ";
1252 
1253  out_stream << '\n' << '\n';
1254 
1255 #endif
1256  }
1257 
1258  }
1259 
1260  // If we wrote any variables, we have to close the variable section now
1261  if (write_variable)
1262  out_stream << "endvars\n";
1263 
1264 
1265  // end of the file
1266  out_stream << "\nendgmv\n";
1267 }
dof_id_type n_active_sub_elem() const
Definition: mesh_base.C:396
double abs(double a)
virtual dof_id_type n_active_elem() const =0
virtual const Point & point(const dof_id_type i) const =0
static UniquePtr< Elem > build(const ElemType type, Elem *p=libmesh_nullptr)
Definition: elem.C:234
unsigned int p_level() const
Definition: elem.h:2218
unsigned int n_partitions() const
Definition: mesh_base.h:811
virtual dof_id_type max_node_id() const =0
The base class for all geometric element types.
Definition: elem.h:86
const class libmesh_nullptr_t libmesh_nullptr
unsigned int & ascii_precision()
IterBase * end
const MT & mesh() const
Definition: mesh_output.h:216
const unsigned int n_vars
Definition: tecplot_io.C:68
bool & partitioning()
Definition: gmv_io.h:107
long double max(long double a, double b)
Base class for Mesh.
Definition: mesh_base.h:67
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
bool & write_subdomain_id_as_material()
Definition: gmv_io.h:114
std::map< std::string, const std::vector< Real > * > _cell_centered_data
Definition: gmv_io.h:225
bool & p_levels()
Definition: gmv_io.h:126
OStreamProxy err(std::cerr)
virtual element_iterator active_elements_begin()=0
std::string enum_to_string(const T e)
virtual element_iterator active_elements_end()=0
bool & subdivide_second_order()
Definition: gmv_io.h:120
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Temporarily serializes a DistributedMesh for output.
virtual unsigned int dim() const =0
virtual unsigned int n_sub_elem() const =0
dof_id_type id() const
Definition: dof_object.h:624
virtual dof_id_type n_nodes() const =0
static ElemType first_order_equivalent_type(const ElemType et)
Definition: elem.C:2746
processor_id_type processor_id()
Definition: libmesh_base.h:96
uint8_t dof_id_type
Definition: id_types.h:64
void libMesh::GMVIO::write_binary ( const std::string &  fname,
const std::vector< Number > *  vec = libmesh_nullptr,
const std::vector< std::string > *  solution_names = libmesh_nullptr 
)
private

This method implements writing a mesh with nodal data to a specified file where the nodal data and variable names are optionally provided.

Definition at line 1275 of file gmv_io.C.

References _cell_centered_data, std::abs(), libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), end, libMesh::err, libMesh::libmesh_assert(), libmesh_nullptr, std::max(), libMesh::MeshInput< MeshBase >::mesh(), libMesh::MeshOutput< MT >::mesh(), libMesh::MeshBase::mesh_dimension(), libMesh::MeshBase::n_active_elem(), libMesh::MeshBase::n_nodes(), libMesh::ParallelObject::n_processors(), n_vars, libMesh::Elem::node_id(), libMesh::Elem::p_level(), p_levels(), partitioning(), libMesh::MeshBase::point(), libMesh::processor_id(), libMesh::DofObject::processor_id(), libMesh::Elem::type(), and write_subdomain_id_as_material().

Referenced by p_levels(), write(), and write_nodal_data().

1278 {
1279  // We use the file-scope global variable eletypes for mapping nodes
1280  // from GMV to libmesh indices, so initialize that data now.
1281  init_eletypes();
1282 
1283  // Get a reference to the mesh
1285 
1286  // This is a parallel_only function
1287  const dof_id_type n_active_elem = mesh.n_active_elem();
1288 
1290  return;
1291 
1292  std::ofstream out_stream (fname.c_str());
1293 
1294  libmesh_assert (out_stream.good());
1295 
1296  unsigned int mesh_max_p_level = 0;
1297 
1298  std::string buffer;
1299 
1300  // Begin interfacing with the GMV data file
1301  {
1302  // write the nodes
1303  buffer = "gmvinput";
1304  out_stream.write(buffer.c_str(), buffer.size());
1305 
1306  buffer = "ieeei4r4";
1307  out_stream.write(buffer.c_str(), buffer.size());
1308  }
1309 
1310 
1311 
1312  // write the nodes
1313  {
1314  buffer = "nodes ";
1315  out_stream.write(buffer.c_str(), buffer.size());
1316 
1317  unsigned int tempint = mesh.n_nodes();
1318  out_stream.write(reinterpret_cast<char *>(&tempint), sizeof(unsigned int));
1319 
1320  // write the x coordinate
1321  std::vector<float> temp(mesh.n_nodes());
1322  for (unsigned int v=0; v<mesh.n_nodes(); v++)
1323  temp[v] = static_cast<float>(mesh.point(v)(0));
1324  out_stream.write(reinterpret_cast<char *>(&temp[0]), sizeof(float)*mesh.n_nodes());
1325 
1326  // write the y coordinate
1327  for (unsigned int v=0; v<mesh.n_nodes(); v++)
1328  {
1329 #if LIBMESH_DIM > 1
1330  temp[v] = static_cast<float>(mesh.point(v)(1));
1331 #else
1332  temp[v] = 0.;
1333 #endif
1334  }
1335  out_stream.write(reinterpret_cast<char *>(&temp[0]), sizeof(float)*mesh.n_nodes());
1336 
1337  // write the z coordinate
1338  for (unsigned int v=0; v<mesh.n_nodes(); v++)
1339  {
1340 #if LIBMESH_DIM > 2
1341  temp[v] = static_cast<float>(mesh.point(v)(2));
1342 #else
1343  temp[v] = 0.;
1344 #endif
1345  }
1346  out_stream.write(reinterpret_cast<char *>(&temp[0]), sizeof(float)*mesh.n_nodes());
1347  }
1348 
1349 
1350  // write the connectivity
1351  {
1352  buffer = "cells ";
1353  out_stream.write(buffer.c_str(), buffer.size());
1354 
1355  unsigned int tempint = n_active_elem;
1356  out_stream.write(reinterpret_cast<char *>(&tempint), sizeof(unsigned int));
1357 
1360 
1361  for ( ; it != end; ++it)
1362  {
1363  const Elem * elem = *it;
1364 
1365  mesh_max_p_level = std::max(mesh_max_p_level,
1366  elem->p_level());
1367 
1368  // The ElementDefinition label contains the GMV name
1369  // and the number of nodes for the respective
1370  // element.
1371  const ElementDefinition & ed = eletypes[elem->type()];
1372 
1373  // The ElementDefinition labels all have a name followed by a
1374  // whitespace and then the number of nodes. To write the
1375  // string in binary, we want to drop everything after that
1376  // space, and then pad the string out to 8 characters.
1377  buffer = ed.label;
1378 
1379  // Erase everything from the first whitespace character to the end of the string.
1380  buffer.erase(buffer.find_first_of(' '), std::string::npos);
1381 
1382  // Append whitespaces until the string is exactly 8 characters long.
1383  while (buffer.size() < 8)
1384  buffer.insert(buffer.end(), ' ');
1385 
1386  // Debugging:
1387  // libMesh::out << "Writing element with name = '" << buffer << "'" << std::endl;
1388 
1389  // Finally, write the 8 character stream to file.
1390  out_stream.write(buffer.c_str(), buffer.size());
1391 
1392  // Write the number of nodes that this type of element has.
1393  // Note: don't want to use the number of nodes of the element,
1394  // because certain elements, like QUAD9 and HEX27 only support
1395  // being written out as lower-order elements (QUAD8 and HEX20,
1396  // respectively).
1397  tempint = ed.node_map.size();
1398  out_stream.write(reinterpret_cast<char *>(&tempint), sizeof(unsigned int));
1399 
1400  // Write the element connectivity
1401  for (std::size_t i=0; i<ed.node_map.size(); i++)
1402  {
1403  dof_id_type id = elem->node_id(ed.node_map[i]) + 1;
1404  out_stream.write(reinterpret_cast<char *>(&id), sizeof(dof_id_type));
1405  }
1406  }
1407  }
1408 
1409 
1410 
1411  // optionally write the partition information
1412  if (this->partitioning())
1413  {
1414  if (this->write_subdomain_id_as_material())
1415  libmesh_error_msg("Not yet supported in GMVIO::write_binary");
1416 
1417  else
1418  {
1419  buffer = "material";
1420  out_stream.write(buffer.c_str(), buffer.size());
1421 
1422  unsigned int tmpint = mesh.n_processors();
1423  out_stream.write(reinterpret_cast<char *>(&tmpint), sizeof(unsigned int));
1424 
1425  tmpint = 0; // IDs are cell based
1426  out_stream.write(reinterpret_cast<char *>(&tmpint), sizeof(unsigned int));
1427 
1428  for (unsigned int proc=0; proc<mesh.n_processors(); proc++)
1429  {
1430  // We have to write exactly 8 bytes. This means that
1431  // the processor id can be up to 3 characters long, and
1432  // we pad it with blank characters on the end.
1433  std::ostringstream oss;
1434  oss << "proc_" << std::setw(3) << std::left << proc;
1435  out_stream.write(oss.str().c_str(), oss.str().size());
1436  }
1437 
1438  std::vector<unsigned int> proc_id (n_active_elem);
1439 
1440  unsigned int n=0;
1441 
1444 
1445  for ( ; it != end; ++it)
1446  {
1447  const Elem * elem = *it;
1448 
1449  // We no longer write sub-elems in GMV, so just assign a
1450  // processor id value to each element.
1451  proc_id[n++] = elem->processor_id() + 1;
1452  }
1453 
1454 
1455  out_stream.write(reinterpret_cast<char *>(&proc_id[0]),
1456  sizeof(unsigned int)*proc_id.size());
1457  }
1458  }
1459 
1460  // If there are *any* variables at all in the system (including
1461  // p level, or arbitrary cell-based data)
1462  // to write, the gmv file needs to contain the word "variable"
1463  // on a line by itself.
1464  bool write_variable = false;
1465 
1466  // 1.) p-levels
1467  if (this->p_levels() && mesh_max_p_level)
1468  write_variable = true;
1469 
1470  // 2.) solution data
1471  if ((solution_names != libmesh_nullptr) && (vec != libmesh_nullptr))
1472  write_variable = true;
1473 
1474  // // 3.) cell-centered data - unsupported
1475  // if (!(this->_cell_centered_data.empty()))
1476  // write_variable = true;
1477 
1478  if (write_variable)
1479  {
1480  buffer = "variable";
1481  out_stream.write(buffer.c_str(), buffer.size());
1482  }
1483 
1484  // optionally write the partition information
1485  if (this->p_levels() && mesh_max_p_level)
1486  {
1487  unsigned int n_floats = n_active_elem;
1488  for (unsigned int i=0; i != mesh.mesh_dimension(); ++i)
1489  n_floats *= 2;
1490 
1491  std::vector<float> temp(n_floats);
1492 
1493  buffer = "p_level";
1494  out_stream.write(buffer.c_str(), buffer.size());
1495 
1496  unsigned int tempint = 0; // p levels are cell data
1497  out_stream.write(reinterpret_cast<char *>(&tempint), sizeof(unsigned int));
1498 
1501  unsigned int n=0;
1502 
1503  for (; it != end; ++it)
1504  for (unsigned int se=0; se<(*it)->n_sub_elem(); se++)
1505  temp[n++] = static_cast<float>( (*it)->p_level() );
1506 
1507  out_stream.write(reinterpret_cast<char *>(&temp[0]),
1508  sizeof(float)*n_floats);
1509  }
1510 
1511 
1512  // optionally write cell-centered data
1513  if (!(this->_cell_centered_data.empty()))
1514  {
1515  libMesh::err << "Cell-centered data not (yet) supported in binary I/O mode!" << std::endl;
1516 
1517  // std::map<std::string, const std::vector<Real> * >::iterator it = this->_cell_centered_data.begin();
1518  // const std::map<std::string, const std::vector<Real> * >::iterator end = this->_cell_centered_data.end();
1519 
1520  // for (; it != end; ++it)
1521  // {
1522  // // Write out the variable name ...
1523  // out_stream.write(it->first.c_str(), it->first.size());
1524 
1525  // // ... followed by a zero.
1526  // unsigned int tempint = 0; // 0 signifies cell data
1527  // out_stream.write(reinterpret_cast<char *>(&tempint), sizeof(unsigned int));
1528 
1529  // // Get a pointer to the array of cell-centered data values
1530  // const std::vector<Real> * the_array = (*it).second;
1531 
1532  // // Since the_array might contain zeros (for inactive elements) we need to
1533  // // make a copy of it containing just values for active elements.
1534  // const unsigned int n_floats = n_active_elem * (1<<mesh.mesh_dimension());
1535  // std::vector<float> temp(n_floats);
1536 
1537  // MeshBase::const_element_iterator elem_it = mesh.active_elements_begin();
1538  // const MeshBase::const_element_iterator elem_end = mesh.active_elements_end();
1539  // unsigned int n=0;
1540 
1541  // for (; elem_it != elem_end; ++elem_it)
1542  // {
1543  // // If there's a seg-fault, it will probably be here!
1544  // const float the_value = static_cast<float>(the_array->operator[]((*elem_it)->id()));
1545 
1546  // for (unsigned int se=0; se<(*elem_it)->n_sub_elem(); se++)
1547  // temp[n++] = the_value;
1548  // }
1549 
1550 
1551  // // Write "the_array" directly to the file
1552  // out_stream.write(reinterpret_cast<char *>(&temp[0]),
1553  // sizeof(float)*n_floats);
1554  // }
1555  }
1556 
1557 
1558 
1559 
1560  // optionally write the data
1561  if ((solution_names != libmesh_nullptr) &&
1562  (vec != libmesh_nullptr))
1563  {
1564  std::vector<float> temp(mesh.n_nodes());
1565 
1566  const unsigned int n_vars =
1567  cast_int<unsigned int>(solution_names->size());
1568 
1569  for (unsigned int c=0; c<n_vars; c++)
1570  {
1571 
1572 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
1573  // for complex data, write three datasets
1574 
1575 
1576  // Real part
1577  buffer = "r_";
1578  out_stream.write(buffer.c_str(), buffer.size());
1579 
1580  buffer = (*solution_names)[c];
1581  out_stream.write(buffer.c_str(), buffer.size());
1582 
1583  unsigned int tempint = 1; // always do nodal data
1584  out_stream.write(reinterpret_cast<char *>(&tempint), sizeof(unsigned int));
1585 
1586  for (unsigned int n=0; n<mesh.n_nodes(); n++)
1587  temp[n] = static_cast<float>( (*vec)[n*n_vars + c].real() );
1588 
1589  out_stream.write(reinterpret_cast<char *>(&temp[0]), sizeof(float)*mesh.n_nodes());
1590 
1591 
1592  // imaginary part
1593  buffer = "i_";
1594  out_stream.write(buffer.c_str(), buffer.size());
1595 
1596  buffer = (*solution_names)[c];
1597  out_stream.write(buffer.c_str(), buffer.size());
1598 
1599  out_stream.write(reinterpret_cast<char *>(&tempint), sizeof(unsigned int));
1600 
1601  for (unsigned int n=0; n<mesh.n_nodes(); n++)
1602  temp[n] = static_cast<float>( (*vec)[n*n_vars + c].imag() );
1603 
1604  out_stream.write(reinterpret_cast<char *>(&temp[0]), sizeof(float)*mesh.n_nodes());
1605 
1606  // magnitude
1607  buffer = "a_";
1608  out_stream.write(buffer.c_str(), buffer.size());
1609  buffer = (*solution_names)[c];
1610  out_stream.write(buffer.c_str(), buffer.size());
1611 
1612  out_stream.write(reinterpret_cast<char *>(&tempint), sizeof(unsigned int));
1613 
1614  for (unsigned int n=0; n<mesh.n_nodes(); n++)
1615  temp[n] = static_cast<float>(std::abs((*vec)[n*n_vars + c]));
1616 
1617  out_stream.write(reinterpret_cast<char *>(&temp[0]), sizeof(float)*mesh.n_nodes());
1618 
1619 #else
1620 
1621  buffer = (*solution_names)[c];
1622  out_stream.write(buffer.c_str(), buffer.size());
1623 
1624  unsigned int tempint = 1; // always do nodal data
1625  out_stream.write(reinterpret_cast<char *>(&tempint), sizeof(unsigned int));
1626 
1627  for (unsigned int n=0; n<mesh.n_nodes(); n++)
1628  temp[n] = static_cast<float>((*vec)[n*n_vars + c]);
1629 
1630  out_stream.write(reinterpret_cast<char *>(&temp[0]), sizeof(float)*mesh.n_nodes());
1631 
1632 #endif
1633  }
1634  }
1635 
1636  // If we wrote any variables, we have to close the variable section now
1637  if (write_variable)
1638  {
1639  buffer = "endvars ";
1640  out_stream.write(buffer.c_str(), buffer.size());
1641  }
1642 
1643  // end the file
1644  buffer = "endgmv ";
1645  out_stream.write(buffer.c_str(), buffer.size());
1646 }
double abs(double a)
virtual dof_id_type n_active_elem() const =0
virtual const Point & point(const dof_id_type i) const =0
virtual ElemType type() const =0
unsigned int p_level() const
Definition: elem.h:2218
processor_id_type n_processors() const
The base class for all geometric element types.
Definition: elem.h:86
const class libmesh_nullptr_t libmesh_nullptr
IterBase * end
const MT & mesh() const
Definition: mesh_output.h:216
const unsigned int n_vars
Definition: tecplot_io.C:68
bool & partitioning()
Definition: gmv_io.h:107
long double max(long double a, double b)
Base class for Mesh.
Definition: mesh_base.h:67
libmesh_assert(j)
bool & write_subdomain_id_as_material()
Definition: gmv_io.h:114
std::map< std::string, const std::vector< Real > * > _cell_centered_data
Definition: gmv_io.h:225
bool & p_levels()
Definition: gmv_io.h:126
OStreamProxy err(std::cerr)
virtual element_iterator active_elements_begin()=0
virtual element_iterator active_elements_end()=0
unsigned int mesh_dimension() const
Definition: mesh_base.C:147
dof_id_type node_id(const unsigned int i) const
Definition: elem.h:1683
virtual dof_id_type n_nodes() const =0
processor_id_type processor_id()
Definition: libmesh_base.h:96
uint8_t dof_id_type
Definition: id_types.h:64
processor_id_type processor_id() const
Definition: dof_object.h:686
void libMesh::GMVIO::write_discontinuous_gmv ( const std::string &  name,
const EquationSystems es,
const bool  write_partitioning,
const std::set< std::string > *  system_names = libmesh_nullptr 
) const

Writes a GMV file with discontinuous data

Definition at line 1656 of file gmv_io.C.

References _cell_centered_data, _write_subdomain_id_as_material, std::abs(), libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), libMesh::EquationSystems::build_discontinuous_solution_vector(), libMesh::EquationSystems::build_variable_names(), libMesh::EDGE2, libMesh::EDGE3, libMesh::EDGE4, end, libMesh::Utility::enum_to_string(), libMesh::err, libMesh::HEX20, libMesh::HEX27, libMesh::HEX8, libMesh::INFEDGE2, libMesh::INFHEX16, libMesh::INFHEX18, libMesh::INFHEX8, libMesh::INFPRISM12, libMesh::INFPRISM6, libMesh::INFQUAD4, libMesh::INFQUAD6, libMesh::libmesh_assert(), libmesh_nullptr, libMesh::MeshInput< MeshBase >::mesh(), libMesh::MeshOutput< MT >::mesh(), libMesh::MeshBase::mesh_dimension(), libMesh::MeshBase::n_active_elem(), libMesh::ParallelObject::n_processors(), n_vars, libMesh::PRISM15, libMesh::PRISM18, libMesh::PRISM6, libMesh::ParallelObject::processor_id(), libMesh::QUAD4, libMesh::QUAD8, libMesh::QUAD9, libMesh::TET10, libMesh::TET4, libMesh::TRI3, and libMesh::TRI6.

Referenced by p_levels(), and libMesh::ErrorVector::plot_error().

1660 {
1661  std::vector<std::string> solution_names;
1662  std::vector<Number> v;
1663 
1664  // Get a reference to the mesh
1666 
1667  es.build_variable_names (solution_names, libmesh_nullptr, system_names);
1668  es.build_discontinuous_solution_vector (v, system_names);
1669 
1670  // These are parallel_only functions
1671  const unsigned int n_active_elem = mesh.n_active_elem();
1672 
1673  if (mesh.processor_id() != 0)
1674  return;
1675 
1676  std::ofstream out_stream(name.c_str());
1677 
1678  libmesh_assert (out_stream.good());
1679 
1680  // Begin interfacing with the GMV data file
1681  {
1682 
1683  // write the nodes
1684  out_stream << "gmvinput ascii" << std::endl << std::endl;
1685 
1686  // Compute the total weight
1687  {
1690 
1691  unsigned int tw=0;
1692 
1693  for ( ; it != end; ++it)
1694  tw += (*it)->n_nodes();
1695 
1696  out_stream << "nodes " << tw << std::endl;
1697  }
1698 
1699 
1700 
1701  // Write all the x values
1702  {
1705 
1706  for ( ; it != end; ++it)
1707  for (unsigned int n=0; n<(*it)->n_nodes(); n++)
1708  out_stream << (*it)->point(n)(0) << " ";
1709 
1710  out_stream << std::endl;
1711  }
1712 
1713 
1714  // Write all the y values
1715  {
1718 
1719  for ( ; it != end; ++it)
1720  for (unsigned int n=0; n<(*it)->n_nodes(); n++)
1721 #if LIBMESH_DIM > 1
1722  out_stream << (*it)->point(n)(1) << " ";
1723 #else
1724  out_stream << 0. << " ";
1725 #endif
1726 
1727  out_stream << std::endl;
1728  }
1729 
1730 
1731  // Write all the z values
1732  {
1735 
1736  for ( ; it != end; ++it)
1737  for (unsigned int n=0; n<(*it)->n_nodes(); n++)
1738 #if LIBMESH_DIM > 2
1739  out_stream << (*it)->point(n)(2) << " ";
1740 #else
1741  out_stream << 0. << " ";
1742 #endif
1743 
1744  out_stream << std::endl << std::endl;
1745  }
1746  }
1747 
1748 
1749 
1750  {
1751  // write the connectivity
1752 
1753  out_stream << "cells " << n_active_elem << std::endl;
1754 
1757 
1758  unsigned int nn=1;
1759 
1760  switch (mesh.mesh_dimension())
1761  {
1762  case 1:
1763  {
1764  for ( ; it != end; ++it)
1765  for (unsigned int se=0; se<(*it)->n_sub_elem(); se++)
1766  {
1767  if (((*it)->type() == EDGE2) ||
1768  ((*it)->type() == EDGE3) ||
1769  ((*it)->type() == EDGE4)
1770 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1771  || ((*it)->type() == INFEDGE2)
1772 #endif
1773  )
1774  {
1775  out_stream << "line 2" << std::endl;
1776  for (unsigned int i=0; i<(*it)->n_nodes(); i++)
1777  out_stream << nn++ << " ";
1778 
1779  }
1780  else
1781  libmesh_error_msg("Unsupported 1D element type: " << Utility::enum_to_string((*it)->type()));
1782 
1783  out_stream << std::endl;
1784  }
1785 
1786  break;
1787  }
1788 
1789  case 2:
1790  {
1791  for ( ; it != end; ++it)
1792  for (unsigned int se=0; se<(*it)->n_sub_elem(); se++)
1793  {
1794  if (((*it)->type() == QUAD4) ||
1795  ((*it)->type() == QUAD8) || // Note: QUAD8 will be output as one central quad and
1796  // four surrounding triangles (though they will be written
1797  // to GMV as QUAD4s).
1798  ((*it)->type() == QUAD9)
1799 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1800  || ((*it)->type() == INFQUAD4)
1801  || ((*it)->type() == INFQUAD6)
1802 #endif
1803  )
1804  {
1805  out_stream << "quad 4" << std::endl;
1806  for (unsigned int i=0; i<(*it)->n_nodes(); i++)
1807  out_stream << nn++ << " ";
1808 
1809  }
1810  else if (((*it)->type() == TRI3) ||
1811  ((*it)->type() == TRI6))
1812  {
1813  out_stream << "tri 3" << std::endl;
1814  for (unsigned int i=0; i<(*it)->n_nodes(); i++)
1815  out_stream << nn++ << " ";
1816 
1817  }
1818  else
1819  libmesh_error_msg("Unsupported 2D element type: " << Utility::enum_to_string((*it)->type()));
1820 
1821  out_stream << std::endl;
1822  }
1823 
1824  break;
1825  }
1826 
1827 
1828  case 3:
1829  {
1830  for ( ; it != end; ++it)
1831  for (unsigned int se=0; se<(*it)->n_sub_elem(); se++)
1832  {
1833  if (((*it)->type() == HEX8) ||
1834  ((*it)->type() == HEX20) ||
1835  ((*it)->type() == HEX27)
1836 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1837  || ((*it)->type() == INFHEX8)
1838  || ((*it)->type() == INFHEX16)
1839  || ((*it)->type() == INFHEX18)
1840 #endif
1841  )
1842  {
1843  out_stream << "phex8 8" << std::endl;
1844  for (unsigned int i=0; i<(*it)->n_nodes(); i++)
1845  out_stream << nn++ << " ";
1846  }
1847  else if (((*it)->type() == PRISM6) ||
1848  ((*it)->type() == PRISM15) ||
1849  ((*it)->type() == PRISM18)
1850 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1851  || ((*it)->type() == INFPRISM6)
1852  || ((*it)->type() == INFPRISM12)
1853 #endif
1854  )
1855  {
1856  out_stream << "pprism6 6" << std::endl;
1857  for (unsigned int i=0; i<(*it)->n_nodes(); i++)
1858  out_stream << nn++ << " ";
1859  }
1860  else if (((*it)->type() == TET4) ||
1861  ((*it)->type() == TET10))
1862  {
1863  out_stream << "tet 4" << std::endl;
1864  for (unsigned int i=0; i<(*it)->n_nodes(); i++)
1865  out_stream << nn++ << " ";
1866  }
1867  else
1868  libmesh_error_msg("Unsupported 3D element type: " << Utility::enum_to_string((*it)->type()));
1869 
1870  out_stream << std::endl;
1871  }
1872 
1873  break;
1874  }
1875 
1876  default:
1877  libmesh_error_msg("Unsupported mesh dimension: " << mesh.mesh_dimension());
1878  }
1879 
1880  out_stream << std::endl;
1881  }
1882 
1883 
1884 
1885  // optionally write the partition information
1886  if (write_partitioning)
1887  {
1889  libmesh_error_msg("Not yet supported in GMVIO::write_discontinuous_gmv");
1890 
1891  else
1892  {
1893  out_stream << "material "
1894  << mesh.n_processors()
1895  << " 0"<< std::endl;
1896 
1897  for (unsigned int proc=0; proc<mesh.n_processors(); proc++)
1898  out_stream << "proc_" << proc << std::endl;
1899 
1902 
1903  for ( ; it != end; ++it)
1904  out_stream << (*it)->processor_id()+1 << std::endl;
1905 
1906  out_stream << std::endl;
1907  }
1908  }
1909 
1910 
1911  // Writing cell-centered data is not yet supported in discontinuous GMV files.
1912  if (!(this->_cell_centered_data.empty()))
1913  {
1914  libMesh::err << "Cell-centered data not (yet) supported for discontinuous GMV files!" << std::endl;
1915  }
1916 
1917 
1918 
1919  // write the data
1920  {
1921  const unsigned int n_vars =
1922  cast_int<unsigned int>(solution_names.size());
1923 
1924  // libmesh_assert_equal_to (v.size(), tw*n_vars);
1925 
1926  out_stream << "variable" << std::endl;
1927 
1928 
1929  for (unsigned int c=0; c<n_vars; c++)
1930  {
1931 
1932 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
1933 
1934  // in case of complex data, write _tree_ data sets
1935  // for each component
1936 
1937  // this is the real part
1938  out_stream << "r_" << solution_names[c] << " 1" << std::endl;
1939  {
1942 
1943  for ( ; it != end; ++it)
1944  for (unsigned int n=0; n<(*it)->n_nodes(); n++)
1945  out_stream << v[(n++)*n_vars + c].real() << " ";
1946  }
1947  out_stream << std::endl << std::endl;
1948 
1949 
1950  // this is the imaginary part
1951  out_stream << "i_" << solution_names[c] << " 1" << std::endl;
1952  {
1955 
1956  for ( ; it != end; ++it)
1957  for (unsigned int n=0; n<(*it)->n_nodes(); n++)
1958  out_stream << v[(n++)*n_vars + c].imag() << " ";
1959  }
1960  out_stream << std::endl << std::endl;
1961 
1962  // this is the magnitude
1963  out_stream << "a_" << solution_names[c] << " 1" << std::endl;
1964  {
1967 
1968  for ( ; it != end; ++it)
1969  for (unsigned int n=0; n<(*it)->n_nodes(); n++)
1970  out_stream << std::abs(v[(n++)*n_vars + c]) << " ";
1971  }
1972  out_stream << std::endl << std::endl;
1973 
1974 #else
1975 
1976  out_stream << solution_names[c] << " 1" << std::endl;
1977  {
1980 
1981  unsigned int nn=0;
1982 
1983  for ( ; it != end; ++it)
1984  for (unsigned int n=0; n<(*it)->n_nodes(); n++)
1985  out_stream << v[(nn++)*n_vars + c] << " ";
1986  }
1987  out_stream << std::endl << std::endl;
1988 
1989 #endif
1990 
1991  }
1992 
1993  out_stream << "endvars" << std::endl;
1994  }
1995 
1996 
1997  // end of the file
1998  out_stream << std::endl << "endgmv" << std::endl;
1999 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
double abs(double a)
virtual dof_id_type n_active_elem() const =0
processor_id_type n_processors() const
const class libmesh_nullptr_t libmesh_nullptr
IterBase * end
const MT & mesh() const
Definition: mesh_output.h:216
const unsigned int n_vars
Definition: tecplot_io.C:68
void build_discontinuous_solution_vector(std::vector< Number > &soln, const std::set< std::string > *system_names=libmesh_nullptr) const
Base class for Mesh.
Definition: mesh_base.h:67
libmesh_assert(j)
std::map< std::string, const std::vector< Real > * > _cell_centered_data
Definition: gmv_io.h:225
bool _write_subdomain_id_as_material
Definition: gmv_io.h:207
OStreamProxy err(std::cerr)
virtual element_iterator active_elements_begin()=0
std::string enum_to_string(const T e)
virtual element_iterator active_elements_end()=0
void build_variable_names(std::vector< std::string > &var_names, const FEType *type=libmesh_nullptr, const std::set< std::string > *system_names=libmesh_nullptr) const
unsigned int mesh_dimension() const
Definition: mesh_base.C:147
processor_id_type processor_id() const
virtual void libMesh::MeshOutput< MeshBase >::write_equation_systems ( const std::string &  ,
const EquationSystems ,
const std::set< std::string > *  system_names = libmesh_nullptr 
)
virtualinherited

This method implements writing a mesh with data to a specified file where the data is taken from the EquationSystems object.

Reimplemented in libMesh::NameBasedIO.

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

void libMesh::GMVIO::write_nodal_data ( const std::string &  fname,
const std::vector< Number > &  soln,
const std::vector< std::string > &  names 
)
virtual

This method implements writing a mesh with nodal data to a specified file where the nodal data and variable names are provided.

Reimplemented from libMesh::MeshOutput< MeshBase >.

Definition at line 277 of file gmv_io.C.

References binary(), write_ascii_old_impl(), and write_binary().

Referenced by libMesh::NameBasedIO::write_nodal_data().

280 {
281  LOG_SCOPE("write_nodal_data()", "GMVIO");
282 
283  if (this->binary())
284  this->write_binary (fname, &soln, &names);
285  else
286  this->write_ascii_old_impl (fname, &soln, &names);
287 }
bool & binary()
Definition: gmv_io.h:95
void write_ascii_old_impl(const std::string &, const std::vector< Number > *=libmesh_nullptr, const std::vector< std::string > *=libmesh_nullptr)
Definition: gmv_io.C:587
void write_binary(const std::string &, const std::vector< Number > *=libmesh_nullptr, const std::vector< std::string > *=libmesh_nullptr)
Definition: gmv_io.C:1275
virtual void libMesh::MeshOutput< MeshBase >::write_nodal_data ( const std::string &  ,
const NumericVector< Number > &  ,
const std::vector< std::string > &   
)
virtualinherited

This method should be overridden by "parallel" output formats for writing nodal data. Instead of getting a localized copy of the nodal solution vector, it is passed a NumericVector of type=PARALLEL which is in node-major order i.e. (u0,v0,w0, u1,v1,w1, u2,v2,w2, u3,v3,w3, ...) and contains n_nodes*n_vars total entries. Then, it is up to the individual I/O class to extract the required solution values from this vector and write them in parallel.

If not implemented, localizes the parallel vector into a std::vector and calls the other version of this function.

Reimplemented in libMesh::Nemesis_IO.

bool& libMesh::GMVIO::write_subdomain_id_as_material ( )
inline

Flag to write element subdomain_id's as GMV "materials" instead of element processor_id's. Allows you to generate exploded views on user-defined subdomains, potentially creating a pretty picture.

Definition at line 114 of file gmv_io.h.

References _write_subdomain_id_as_material.

Referenced by write_ascii_new_impl(), write_ascii_old_impl(), and write_binary().

bool _write_subdomain_id_as_material
Definition: gmv_io.h:207

Member Data Documentation

bool libMesh::GMVIO::_binary
private

Flag to write binary data.

Definition at line 191 of file gmv_io.h.

Referenced by binary().

std::map<std::string, const std::vector<Real> * > libMesh::GMVIO::_cell_centered_data
private

Storage for arbitrary cell-centered data. Ex: You can use this to plot the effectivity index for a given cell. The map is between the string representing the variable name and a pointer to a vector containing the data.

Definition at line 225 of file gmv_io.h.

Referenced by add_cell_centered_data(), write_ascii_new_impl(), write_ascii_old_impl(), write_binary(), and write_discontinuous_gmv().

bool libMesh::GMVIO::_discontinuous
private

Flag to write the mesh as discontinuous patches.

Definition at line 196 of file gmv_io.h.

Referenced by discontinuous().

const bool libMesh::MeshOutput< MeshBase >::_is_parallel_format
protectedinherited

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 141 of file mesh_output.h.

Referenced by libMesh::FroIO::write(), libMesh::PostscriptIO::write(), and libMesh::EnsightIO::write().

unsigned int libMesh::GMVIO::_next_elem_id
private

Definition at line 231 of file gmv_io.h.

Referenced by _read_one_cell(), and read().

std::map<std::string, std::vector<Number> > libMesh::GMVIO::_nodal_data
private

Definition at line 236 of file gmv_io.h.

Referenced by _read_var(), and copy_nodal_solution().

bool libMesh::GMVIO::_p_levels
private

Flag to write the mesh p refinement levels.

Definition at line 217 of file gmv_io.h.

Referenced by p_levels().

bool libMesh::GMVIO::_partitioning
private

Flag to write the mesh partitioning.

Definition at line 201 of file gmv_io.h.

Referenced by partitioning().

std::map< std::string, ElemType > libMesh::GMVIO::_reading_element_map = GMVIO::build_reading_element_map()
staticprivate

Static map from string -> ElementType for use during reading.

Definition at line 242 of file gmv_io.h.

Referenced by gmv_elem_to_libmesh_elem().

const bool libMesh::MeshOutput< MeshBase >::_serial_only_needed_on_proc_0
protectedinherited

Flag specifying whether this format can be written by only serializing the mesh to processor zero

If this is false (default) the mesh will be serialized to all processors

Definition at line 150 of file mesh_output.h.

bool libMesh::GMVIO::_subdivide_second_order
private

Flag to subdivide second order elements

Definition at line 212 of file gmv_io.h.

Referenced by subdivide_second_order().

bool libMesh::GMVIO::_write_subdomain_id_as_material
private

Flag to write element subdomain_id's as GMV "materials" instead of element processor_id's.

Definition at line 207 of file gmv_io.h.

Referenced by write_discontinuous_gmv(), and write_subdomain_id_as_material().


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