libMesh
Public Member Functions | Protected Member Functions | Protected Attributes | Private Member Functions | Static Private Member Functions | Private Attributes | Static Private Attributes | List of all members
libMesh::GMVIO Class Reference

This class implements writing meshes in the GMV format. More...

#include <gmv_io.h>

Inheritance diagram for libMesh::GMVIO:
[legend]

Public Member Functions

 GMVIO (const MeshBase &)
 Constructor. More...
 
 GMVIO (MeshBase &)
 Constructor. More...
 
virtual void write (const std::string &) libmesh_override
 This method implements writing a mesh to a specified file. More...
 
virtual void read (const std::string &mesh_file) libmesh_override
 This method implements reading a mesh from a specified file. More...
 
virtual void write_nodal_data (const std::string &, const std::vector< Number > &, const std::vector< std::string > &) libmesh_override
 This method implements writing a mesh with nodal data to a specified file where the nodal data and variable names are provided. More...
 
bool & binary ()
 Flag indicating whether or not to write a binary file. More...
 
bool & discontinuous ()
 Flag indicating whether or not to write the mesh as discontinuous cell patches. More...
 
bool & partitioning ()
 Flag indicating whether or not to write the partitioning information for the mesh. More...
 
bool & write_subdomain_id_as_material ()
 Flag to write element subdomain_id's as GMV "materials" instead of element processor_id's. More...
 
bool & subdivide_second_order ()
 Flag indicating whether or not to subdivide second order elements. More...
 
bool & p_levels ()
 Flag indicating whether or not to write p level information for p refined meshes. More...
 
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
 Writes a GMV file with discontinuous data. More...
 
void write_ascii_new_impl (const std::string &, const std::vector< Number > *=libmesh_nullptr, const std::vector< std::string > *=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. More...
 
void 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. More...
 
void 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. More...
 
virtual void write_equation_systems (const std::string &, const EquationSystems &, const std::set< std::string > *system_names=libmesh_nullptr)
 This method implements writing a mesh with data to a specified file where the data is taken from the EquationSystems object. More...
 
virtual void write_nodal_data (const std::string &, const NumericVector< Number > &, const std::vector< std::string > &)
 This method should be overridden by "parallel" output formats for writing nodal data. More...
 
unsigned intascii_precision ()
 Return/set the precision to use when writing ASCII files. More...
 

Protected Member Functions

MeshBasemesh ()
 
void set_n_partitions (unsigned int n_parts)
 Sets the number of partitions in the mesh. More...
 
void skip_comment_lines (std::istream &in, const char comment_start)
 Reads input from in, skipping all the lines that start with the character comment_start. More...
 
const MeshBasemesh () const
 

Protected Attributes

std::vector< bool > elems_of_dimension
 A vector of bools describing what dimension elements have been encountered when reading a mesh. More...
 
const bool _is_parallel_format
 Flag specifying whether this format is parallel-capable. More...
 
const bool _serial_only_needed_on_proc_0
 Flag specifying whether this format can be written by only serializing the mesh to processor zero. More...
 

Private Member Functions

void write_ascii_old_impl (const std::string &, const std::vector< Number > *=libmesh_nullptr, const std::vector< std::string > *=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. More...
 
void write_binary (const std::string &, const std::vector< Number > *=libmesh_nullptr, const std::vector< std::string > *=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. More...
 
void _read_nodes ()
 Helper functions for reading nodes/cells from a GMV file. More...
 
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 ()
 Static function used to build the _reading_element_map. More...
 

Private Attributes

bool _binary
 Flag to write binary data. More...
 
bool _discontinuous
 Flag to write the mesh as discontinuous patches. More...
 
bool _partitioning
 Flag to write the mesh partitioning. More...
 
bool _write_subdomain_id_as_material
 Flag to write element subdomain_id's as GMV "materials" instead of element processor_id's. More...
 
bool _subdivide_second_order
 Flag to subdivide second order elements. More...
 
bool _p_levels
 Flag to write the mesh p refinement levels. More...
 
std::map< std::string, const std::vector< Real > * > _cell_centered_data
 Storage for arbitrary cell-centered data. More...
 
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()
 Static map from string -> ElementType for use during reading. More...
 

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
Flag to write binary data.
Definition: gmv_io.h:191
unsigned int _next_elem_id
Definition: gmv_io.h:231
bool _discontinuous
Flag to write the mesh as discontinuous patches.
Definition: gmv_io.h:196
bool _write_subdomain_id_as_material
Flag to write element subdomain_id&#39;s as GMV "materials" instead of element processor_id&#39;s.
Definition: gmv_io.h:207
bool _subdivide_second_order
Flag to subdivide second order elements.
Definition: gmv_io.h:212
bool _p_levels
Flag to write the mesh p refinement levels.
Definition: gmv_io.h:217
bool _partitioning
Flag to write the mesh 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
Flag to write binary data.
Definition: gmv_io.h:191
unsigned int _next_elem_id
Definition: gmv_io.h:231
bool _discontinuous
Flag to write the mesh as discontinuous patches.
Definition: gmv_io.h:196
bool _write_subdomain_id_as_material
Flag to write element subdomain_id&#39;s as GMV "materials" instead of element processor_id&#39;s.
Definition: gmv_io.h:207
bool _subdivide_second_order
Flag to subdivide second order elements.
Definition: gmv_io.h:212
bool _p_levels
Flag to write the mesh p refinement levels.
Definition: gmv_io.h:217
bool _partitioning
Flag to write the mesh partitioning.
Definition: gmv_io.h:201

Member Function Documentation

void libMesh::GMVIO::_read_materials ( )
private

Definition at line 2055 of file gmv_io.C.

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

Referenced by read().

2056 {
2057 #ifdef LIBMESH_HAVE_GMV
2058 
2059  // LibMesh assigns materials on a per-cell basis
2060  libmesh_assert_equal_to (GMVLib::gmv_data.datatype, CELL);
2061 
2062  // // Material names: LibMesh has no use for these currently...
2063  // libMesh::out << "Number of material names="
2064  // << GMVLib::gmv_data.num
2065  // << std::endl;
2066 
2067  // for (int i = 0; i < GMVLib::gmv_data.num; i++)
2068  // {
2069  // // Build a 32-char string from the appropriate entries
2070  // std::string mat_string(&GMVLib::gmv_data.chardata1[i*33], 32);
2071 
2072  // libMesh::out << "Material name " << i << ": " << mat_string << std::endl;
2073  // }
2074 
2075  // // Material labels: These correspond to (1-based) CPU IDs, and
2076  // // there should be 1 of these for each element.
2077  // libMesh::out << "Number of material labels = "
2078  // << GMVLib::gmv_data.nlongdata1
2079  // << std::endl;
2080 
2081  for (int i = 0; i < GMVLib::gmv_data.nlongdata1; i++)
2082  {
2083  // Debugging Info
2084  // libMesh::out << "Material ID " << i << ": "
2085  // << GMVLib::gmv_data.longdata1[i]
2086  // << std::endl;
2087 
2088  MeshInput<MeshBase>::mesh().elem_ref(i).processor_id() =
2089  cast_int<processor_id_type>(GMVLib::gmv_data.longdata1[i]-1);
2090  }
2091 
2092 #endif
2093 }
void libMesh::GMVIO::_read_nodes ( )
private

Helper functions for reading nodes/cells from a GMV file.

Definition at line 2098 of file gmv_io.C.

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

Referenced by read().

2099 {
2100 #ifdef LIBMESH_HAVE_GMV
2101  // Debugging
2102  // libMesh::out << "gmv_data.datatype = " << GMVLib::gmv_data.datatype << std::endl;
2103 
2104  // LibMesh writes UNSTRUCT=100 node data
2105  libmesh_assert_equal_to (GMVLib::gmv_data.datatype, UNSTRUCT);
2106 
2107  // The nodal data is stored in gmv_data.doubledata{1,2,3}
2108  // and is nnodes long
2109  for (int i = 0; i < GMVLib::gmv_data.num; i++)
2110  {
2111  // Debugging
2112  // libMesh::out << "(x,y,z)="
2113  // << "("
2114  // << GMVLib::gmv_data.doubledata1[i] << ","
2115  // << GMVLib::gmv_data.doubledata2[i] << ","
2116  // << GMVLib::gmv_data.doubledata3[i]
2117  // << ")"
2118  // << std::endl;
2119 
2120  // Add the point to the Mesh
2121  MeshInput<MeshBase>::mesh().add_point(Point(GMVLib::gmv_data.doubledata1[i],
2122  GMVLib::gmv_data.doubledata2[i],
2123  GMVLib::gmv_data.doubledata3[i]), i);
2124  }
2125 #endif
2126 }
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
void libMesh::GMVIO::_read_one_cell ( )
private

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

2130 {
2131 #ifdef LIBMESH_HAVE_GMV
2132  // Debugging
2133  // libMesh::out << "gmv_data.datatype=" << GMVLib::gmv_data.datatype << std::endl;
2134 
2135  // This is either a REGULAR=111 cell or
2136  // the ENDKEYWORD=207 of the cells
2137 #ifndef NDEBUG
2138  bool recognized =
2139  (GMVLib::gmv_data.datatype==REGULAR) ||
2140  (GMVLib::gmv_data.datatype==ENDKEYWORD);
2141 #endif
2142  libmesh_assert (recognized);
2143 
2145 
2146  if (GMVLib::gmv_data.datatype == REGULAR)
2147  {
2148  // Debugging
2149  // libMesh::out << "Name of the cell is: " << GMVLib::gmv_data.name1 << std::endl;
2150  // libMesh::out << "Cell has " << GMVLib::gmv_data.num2 << " vertices." << std::endl;
2151 
2152  // We need a mapping from GMV element types to LibMesh
2153  // ElemTypes. Basically the reverse of the eletypes
2154  // std::map above.
2155  //
2156  // FIXME: Since Quad9's apparently don't exist for GMV, and since
2157  // In general we write linear sub-elements to GMV files, we need
2158  // to be careful to read back in exactly what we wrote out...
2159  //
2160  // The gmv_data.name1 field is padded with blank characters when
2161  // it is read in by GMV, so the function that converts it to the
2162  // libmesh element type needs to take that into account.
2163  ElemType type = this->gmv_elem_to_libmesh_elem(GMVLib::gmv_data.name1);
2164 
2165  Elem * elem = Elem::build(type).release();
2166  elem->set_id(_next_elem_id++);
2167 
2168  // Get the ElementDefinition object for this element type
2169  const ElementDefinition & eledef = eletypes[type];
2170 
2171  // Print out the connectivity information for
2172  // this cell.
2173  for (int i=0; i<GMVLib::gmv_data.num2; i++)
2174  {
2175  // Debugging
2176  // libMesh::out << "Vertex " << i << " is node " << GMVLib::gmv_data.longdata1[i] << std::endl;
2177 
2178  // Map index i to GMV's numbering scheme
2179  unsigned mapped_i = eledef.node_map[i];
2180 
2181  // Note: Node numbers (as stored in libmesh) are 1-based
2182  elem->set_node(i) = mesh.node_ptr
2183  (cast_int<dof_id_type>(GMVLib::gmv_data.longdata1[mapped_i]-1));
2184  }
2185 
2186  elems_of_dimension[elem->dim()] = true;
2187 
2188  // Add the newly-created element to the mesh
2189  mesh.add_elem(elem);
2190  }
2191 
2192 
2193  if (GMVLib::gmv_data.datatype == ENDKEYWORD)
2194  {
2195  // There isn't a cell to read, so we just return
2196  return;
2197  }
2198 
2199 #endif
2200 }
virtual Node *& set_node(const unsigned int i)
Definition: elem.h:1941
static UniquePtr< Elem > build(const ElemType type, Elem *p=libmesh_nullptr)
Definition: elem.C:238
std::vector< bool > elems_of_dimension
A vector of bools describing what dimension elements have been encountered when reading a mesh...
Definition: mesh_input.h:97
ElemType
Defines an enum for geometric element types.
unsigned int _next_elem_id
Definition: gmv_io.h:231
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
virtual const Node * node_ptr(const dof_id_type i) const =0
This is the MeshBase class.
Definition: mesh_base.h:68
dof_id_type & set_id()
Definition: dof_object.h:641
libmesh_assert(j)
ElemType gmv_elem_to_libmesh_elem(std::string elemname)
Definition: gmv_io.C:2203
virtual Elem * add_elem(Elem *e)=0
Add elem e to the end of the element array.
virtual unsigned int dim() const =0
void libMesh::GMVIO::_read_var ( )
private

Definition at line 2043 of file gmv_io.C.

References _nodal_data.

Referenced by read().

2044 {
2045 #ifdef LIBMESH_HAVE_GMV
2046 
2047  // Copy all the variable's values into a local storage vector.
2048  _nodal_data.insert ( std::make_pair(std::string(GMVLib::gmv_data.name1),
2049  std::vector<Number>(GMVLib::gmv_data.doubledata1, GMVLib::gmv_data.doubledata1+GMVLib::gmv_data.num) ) );
2050 #endif
2051 }
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 1869 of file gmv_io.C.

References _cell_centered_data, and libMesh::libmesh_assert().

Referenced by p_levels().

1871 {
1872  libmesh_assert(cell_centered_data_vals);
1873 
1874  // Make sure there are *at least* enough entries for all the active elements.
1875  // There can also be entries for inactive elements, they will be ignored.
1876  // libmesh_assert_greater_equal (cell_centered_data_vals->size(),
1877  // MeshOutput<MeshBase>::mesh().n_active_elem());
1878  this->_cell_centered_data[cell_centered_data_name] = cell_centered_data_vals;
1879 }
libmesh_assert(j)
std::map< std::string, const std::vector< Real > * > _cell_centered_data
Storage for arbitrary 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 ( )

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
Flag to write binary data.
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 2220 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().

2221 {
2222  // Check for easy return if there isn't any nodal data
2223  if (_nodal_data.empty())
2224  {
2225  libMesh::err << "Unable to copy nodal solution: No nodal "
2226  << "solution has been read in from file." << std::endl;
2227  return;
2228  }
2229 
2230  // Be sure there is at least one system
2231  libmesh_assert (es.n_systems());
2232 
2233  // Keep track of variable names which have been found and
2234  // copied already. This could be used to prevent us from
2235  // e.g. copying the same var into 2 different systems ...
2236  // but this seems unlikely. Also, it is used to tell if
2237  // any variables which were read in were not successfully
2238  // copied to the EquationSystems.
2239  std::set<std::string> vars_copied;
2240 
2241  // For each entry in the nodal data map, try to find a system
2242  // that has the same variable key name.
2243  for (unsigned int sys=0; sys<es.n_systems(); ++sys)
2244  {
2245  // Get a generic reference to the current System
2246  System & system = es.get_system(sys);
2247 
2248  // And a reference to that system's dof_map
2249  // const DofMap & dof_map = system.get_dof_map();
2250 
2251  // For each var entry in the _nodal_data map, try to find
2252  // that var in the system
2253  std::map<std::string, std::vector<Number>>::iterator it = _nodal_data.begin();
2254  const std::map<std::string, std::vector<Number>>::iterator end = _nodal_data.end();
2255  for (; it != end; ++it)
2256  {
2257  std::string var_name = it->first;
2258  // libMesh::out << "Searching for var " << var_name << " in system " << sys << std::endl;
2259 
2260  if (system.has_variable(var_name))
2261  {
2262  // Check if there are as many nodes in the mesh as there are entries
2263  // in the stored nodal data vector
2264  libmesh_assert_equal_to ( it->second.size(), MeshInput<MeshBase>::mesh().n_nodes() );
2265 
2266  const unsigned int var_num = system.variable_number(var_name);
2267 
2268  // libMesh::out << "Variable "
2269  // << var_name
2270  // << " is variable "
2271  // << var_num
2272  // << " in system " << sys << std::endl;
2273 
2274  // The only type of nodal data we can read in from GMV is for
2275  // linear LAGRANGE type elements.
2276  const FEType & fe_type = system.variable_type(var_num);
2277  if ((fe_type.order.get_order() != FIRST) || (fe_type.family != LAGRANGE))
2278  {
2279  libMesh::err << "Only FIRST-order LAGRANGE variables can be read from GMV files. "
2280  << "Skipping variable " << var_name << std::endl;
2281  break;
2282  }
2283 
2284 
2285  // Loop over the stored vector's entries, inserting them into
2286  // the System's solution if appropriate.
2287  for (std::size_t i=0; i<it->second.size(); ++i)
2288  {
2289  // Since this var came from a GMV file, the index i corresponds to
2290  // the (single) DOF value of the current variable for node i.
2291  const unsigned int dof_index =
2292  MeshInput<MeshBase>::mesh().node_ptr(i)->dof_number(sys, // system #
2293  var_num, // var #
2294  0); // component #, always zero for LAGRANGE
2295 
2296  // libMesh::out << "Value " << i << ": "
2297  // << it->second [i]
2298  // << ", dof index="
2299  // << dof_index << std::endl;
2300 
2301  // If the dof_index is local to this processor, set the value
2302  if ((dof_index >= system.solution->first_local_index()) &&
2303  (dof_index < system.solution->last_local_index()))
2304  system.solution->set (dof_index, it->second [i]);
2305  } // end loop over my GMVIO's copy of the solution
2306 
2307  // Add the most recently copied var to the set of copied vars
2308  vars_copied.insert (var_name);
2309  } // end if (system.has_variable)
2310  } // end for loop over _nodal_data
2311 
2312  // Communicate parallel values before going to the next system.
2313  system.solution->close();
2314  system.update();
2315 
2316  } // end loop over all systems
2317 
2318 
2319 
2320  // Warn the user if any GMV variables were not successfully copied over to the EquationSystems object
2321  {
2322  std::map<std::string, std::vector<Number>>::iterator it = _nodal_data.begin();
2323  const std::map<std::string, std::vector<Number>>::iterator end = _nodal_data.end();
2324 
2325  for (; it != end; ++it)
2326  {
2327  if (vars_copied.find( it->first ) == vars_copied.end())
2328  {
2329  libMesh::err << "Warning: Variable "
2330  << it->first
2331  << " was not copied to the EquationSystems object."
2332  << std::endl;
2333  }
2334  }
2335  }
2336 
2337 }
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:178
OStreamProxy err
FEFamily family
The type of finite element.
Definition: fe_type.h:203
int get_order() const
Explicitly request the order as an int.
Definition: fe_type.h:77
ImplicitSystem & sys
bool has_variable(const std::string &var) const
Definition: system.C:1256
OrderWrapper order
The approximation order of the element.
Definition: fe_type.h:197
IterBase * end
Also have a polymorphic pointer to the end object, this prevents iterating past the end...
const FEType & variable_type(const unsigned int i) const
Definition: system.h:2164
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:1263
This is the base class for classes which contain information related to any physical process that mig...
Definition: system.h:76
UniquePtr< NumericVector< Number > > solution
Data structure to hold solution values.
Definition: system.h:1523
unsigned int n_systems() const
virtual void update()
Update the local values to reflect the solution on neighboring processors.
Definition: system.C:425
const T_sys & get_system(const std::string &name) const
bool& libMesh::GMVIO::discontinuous ( )

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
Flag to write the mesh as discontinuous patches.
Definition: gmv_io.h:196
ElemType libMesh::GMVIO::gmv_elem_to_libmesh_elem ( std::string  elemname)
private

Definition at line 2203 of file gmv_io.C.

References _reading_element_map.

Referenced by _read_one_cell().

2204 {
2205  // Erase any whitespace padding in name coming from gmv before performing comparison.
2206  elemname.erase(std::remove_if(elemname.begin(), elemname.end(), isspace), elemname.end());
2207 
2208  // Look up the string in our string->ElemType name.
2209  std::map<std::string, ElemType>::iterator it = _reading_element_map.find(elemname);
2210 
2211  if (it == _reading_element_map.end())
2212  libmesh_error_msg("Unknown/unsupported element: " << elemname << " was read.");
2213 
2214  return it->second;
2215 }
static std::map< std::string, ElemType > _reading_element_map
Static map from string -> ElementType for use during reading.
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::CheckpointIO::read_header(), libMesh::XdrIO::read_header(), 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().

const MeshBase & libMesh::MeshOutput< MeshBase >::mesh ( ) const
protectedinherited
bool& libMesh::GMVIO::p_levels ( )

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
Flag to write the mesh p refinement levels.
Definition: gmv_io.h:217
bool& libMesh::GMVIO::partitioning ( )

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
Flag to write the mesh 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 1886 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().

1887 {
1888  // This is a serial-only process for now;
1889  // the Mesh should be read on processor 0 and
1890  // broadcast later
1891  libmesh_assert_equal_to (MeshOutput<MeshBase>::mesh().processor_id(), 0);
1892 
1893  _next_elem_id = 0;
1894 
1895  libmesh_experimental();
1896 
1897 #ifndef LIBMESH_HAVE_GMV
1898 
1899  libmesh_error_msg("Cannot read GMV file " << name << " without the GMV API.");
1900 
1901 #else
1902  // We use the file-scope global variable eletypes for mapping nodes
1903  // from GMV to libmesh indices, so initialize that data now.
1904  init_eletypes();
1905 
1906  // Clear the mesh so we are sure to start from a pristine state.
1908  mesh.clear();
1909 
1910  // Keep track of what kinds of elements this file contains
1911  elems_of_dimension.clear();
1912  elems_of_dimension.resize(4, false);
1913 
1914  // It is apparently possible for gmv files to contain
1915  // a "fromfile" directive (?) But we currently don't make
1916  // any use of this feature in LibMesh. Nonzero return val
1917  // from any function usually means an error has occurred.
1918  int ierr = GMVLib::gmvread_open_fromfileskip(const_cast<char *>(name.c_str()));
1919  if (ierr != 0)
1920  libmesh_error_msg("GMVLib::gmvread_open_fromfileskip failed!");
1921 
1922 
1923  // Loop through file until GMVEND.
1924  int iend = 0;
1925  while (iend == 0)
1926  {
1927  GMVLib::gmvread_data();
1928 
1929  // Check for GMVEND.
1930  if (GMVLib::gmv_data.keyword == GMVEND)
1931  {
1932  iend = 1;
1933  GMVLib::gmvread_close();
1934  break;
1935  }
1936 
1937  // Check for GMVERROR.
1938  if (GMVLib::gmv_data.keyword == GMVERROR)
1939  libmesh_error_msg("Encountered GMVERROR while reading!");
1940 
1941  // Process the data.
1942  switch (GMVLib::gmv_data.keyword)
1943  {
1944  case NODES:
1945  {
1946  //libMesh::out << "Reading nodes." << std::endl;
1947 
1948  if (GMVLib::gmv_data.num2 == NODES)
1949  this->_read_nodes();
1950 
1951  else if (GMVLib::gmv_data.num2 == NODE_V)
1952  libmesh_error_msg("Unsupported GMV data type NODE_V!");
1953 
1954  break;
1955  }
1956 
1957  case CELLS:
1958  {
1959  // Read 1 cell at a time
1960  // libMesh::out << "\nReading one cell." << std::endl;
1961  this->_read_one_cell();
1962  break;
1963  }
1964 
1965  case MATERIAL:
1966  {
1967  // keyword == 6
1968  // These are the materials, which we use to specify the mesh
1969  // partitioning.
1970  this->_read_materials();
1971  break;
1972  }
1973 
1974  case VARIABLE:
1975  {
1976  // keyword == 8
1977  // This is a field variable.
1978 
1979  // Check to see if we're done reading variables and break out.
1980  if (GMVLib::gmv_data.datatype == ENDKEYWORD)
1981  {
1982  // libMesh::out << "Done reading GMV variables." << std::endl;
1983  break;
1984  }
1985 
1986  if (GMVLib::gmv_data.datatype == NODE)
1987  {
1988  // libMesh::out << "Reading node field data for variable "
1989  // << GMVLib::gmv_data.name1 << std::endl;
1990  this->_read_var();
1991  break;
1992  }
1993 
1994  else
1995  {
1996  libMesh::err << "Warning: Skipping variable: "
1997  << GMVLib::gmv_data.name1
1998  << " which is of unsupported GMV datatype "
1999  << GMVLib::gmv_data.datatype
2000  << ". Nodal field data is currently the only type currently supported."
2001  << std::endl;
2002  break;
2003  }
2004 
2005  }
2006 
2007  default:
2008  libmesh_error_msg("Encountered unknown GMV keyword " << GMVLib::gmv_data.keyword);
2009 
2010  } // end switch
2011  } // end while
2012 
2013  // Set the mesh dimension to the largest encountered for an element
2014  for (unsigned char i=0; i!=4; ++i)
2015  if (elems_of_dimension[i])
2016  mesh.set_mesh_dimension(i);
2017 
2018 #if LIBMESH_DIM < 3
2019  if (mesh.mesh_dimension() > LIBMESH_DIM)
2020  libmesh_error_msg("Cannot open dimension " \
2021  << mesh.mesh_dimension() \
2022  << " mesh file when configured without " \
2023  << mesh.mesh_dimension() \
2024  << "D support.");
2025 #endif
2026 
2027  // Done reading in the mesh, now call find_neighbors, etc.
2028  // mesh.find_neighbors();
2029 
2030  // Don't allow renumbering of nodes and elements when calling
2031  // prepare_for_use(). It is usually confusing for users when the
2032  // Mesh gets renumbered automatically, since sometimes there are
2033  // other parts of the simulation which rely on the original
2034  // ordering.
2035  mesh.allow_renumbering(false);
2036  mesh.prepare_for_use();
2037 #endif
2038 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:39
OStreamProxy err
void _read_var()
Definition: gmv_io.C:2043
void _read_nodes()
Helper functions for reading nodes/cells from a GMV file.
Definition: gmv_io.C:2098
void allow_renumbering(bool allow)
If false is passed in then this mesh will no longer be renumbered when being prepared for use...
Definition: mesh_base.h:749
std::vector< bool > elems_of_dimension
A vector of bools describing what dimension elements have been encountered when reading a mesh...
Definition: mesh_input.h:97
unsigned int _next_elem_id
Definition: gmv_io.h:231
This class defines an abstract interface for Mesh output.
Definition: mesh_output.h:53
This is the MeshBase class.
Definition: mesh_base.h:68
void prepare_for_use(const bool skip_renumber_nodes_and_elements=false, const bool skip_find_neighbors=false)
Prepare a newly created (or read) mesh for use.
Definition: mesh_base.C:174
void set_mesh_dimension(unsigned char d)
Resets the logical dimension of the mesh.
Definition: mesh_base.h:199
void _read_materials()
Definition: gmv_io.C:2055
virtual void clear()
Deletes all the data that are currently stored.
Definition: mesh_base.C:285
PetscErrorCode ierr
unsigned int mesh_dimension() const
Definition: mesh_base.C:148
void _read_one_cell()
Definition: gmv_io.C:2129
processor_id_type processor_id()
Definition: libmesh_base.h:96
void libMesh::MeshInput< MeshBase >::set_n_partitions ( unsigned int  n_parts)
protectedinherited

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

91 { this->mesh().set_n_partitions() = n_parts; }
unsigned int & set_n_partitions()
Definition: mesh_base.h:1329
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 ( )

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
Flag to subdivide second order elements.
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 main(), and 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()
Flag indicating whether or not to write a binary file.
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)
This method implements writing a mesh with nodal data to a specified file where the nodal data and va...
Definition: gmv_io.C:569
void write_binary(const std::string &, const std::vector< Number > *=libmesh_nullptr, const std::vector< std::string > *=libmesh_nullptr)
This method implements writing a mesh with nodal data to a specified file where the nodal data and va...
Definition: gmv_io.C:1232
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_element_ptr_range(), libMesh::MeshOutput< MeshBase >::ascii_precision(), end, libMesh::err, 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::MeshBase::n_partitions(), n_vars, p_levels(), partitioning(), libMesh::MeshBase::point(), libMesh::processor_id(), libMesh::Real, subdivide_second_order(), 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 
360  for (const auto & elem : mesh.active_element_ptr_range())
361  {
362  mesh_max_p_level = std::max(mesh_max_p_level,
363  elem->p_level());
364 
365  // Make sure we have a valid entry for
366  // the current element type.
367  libmesh_assert (eletypes.count(elem->type()));
368 
369  const ElementDefinition & ele = eletypes[elem->type()];
370 
371  // The element mapper better not require any more nodes
372  // than are present in the current element!
373  libmesh_assert_less_equal (ele.node_map.size(), elem->n_nodes());
374 
375  out_stream << ele.label << "\n";
376  for (std::size_t i=0; i < ele.node_map.size(); i++)
377  out_stream << elem->node_id(ele.node_map[i])+1 << " ";
378  out_stream << "\n";
379  }
380  out_stream << "\n";
381  }
382 
383  // optionally write the partition information
384  if (this->partitioning())
385  {
386  if (this->write_subdomain_id_as_material())
387  libmesh_error_msg("Not yet supported in GMVIO::write_ascii_new_impl");
388 
389  else // write processor IDs as materials. This is the default
390  {
391  out_stream << "material "
392  << mesh.n_partitions()
393  // Note: GMV may give you errors like
394  // Error, material for cell 1 is greater than 1
395  // Error, material for cell 2 is greater than 1
396  // Error, material for cell 3 is greater than 1
397  // ... because you put the wrong number of partitions here.
398  // To ensure you write the correct number of materials, call
399  // mesh.recalculate_n_partitions() before writing out the
400  // mesh.
401  // Note: we can't call it now because the Mesh is const here and
402  // it is a non-const function.
403  << " 0\n";
404 
405  for (unsigned int proc=0; proc<mesh.n_partitions(); proc++)
406  out_stream << "proc_" << proc << "\n";
407 
408  // FIXME - don't we need to use an ElementDefinition here? - RHS
409  for (const auto & elem : mesh.active_element_ptr_range())
410  out_stream << elem->processor_id()+1 << "\n";
411  out_stream << "\n";
412  }
413  }
414 
415  // If there are *any* variables at all in the system (including
416  // p level, or arbitrary cell-based data)
417  // to write, the gmv file needs to contain the word "variable"
418  // on a line by itself.
419  bool write_variable = false;
420 
421  // 1.) p-levels
422  if (this->p_levels() && mesh_max_p_level)
423  write_variable = true;
424 
425  // 2.) solution data
426  if ((solution_names != libmesh_nullptr) && (v != libmesh_nullptr))
427  write_variable = true;
428 
429  // 3.) cell-centered data
430  if (!(this->_cell_centered_data.empty()))
431  write_variable = true;
432 
433  if (write_variable)
434  out_stream << "variable\n";
435 
436  // if ((this->p_levels() && mesh_max_p_level) ||
437  // ((solution_names != libmesh_nullptr) && (v != libmesh_nullptr)))
438  // out_stream << "variable\n";
439 
440  // optionally write the polynomial degree information
441  if (this->p_levels() && mesh_max_p_level)
442  {
443  out_stream << "p_level 0\n";
444 
445  for (const auto & elem : mesh.active_element_ptr_range())
446  {
447  const ElementDefinition & ele = eletypes[elem->type()];
448 
449  // The element mapper better not require any more nodes
450  // than are present in the current element!
451  libmesh_assert_less_equal (ele.node_map.size(), elem->n_nodes());
452 
453  for (std::size_t i=0; i < ele.node_map.size(); i++)
454  out_stream << elem->p_level() << " ";
455  }
456  out_stream << "\n\n";
457  }
458 
459 
460  // optionally write cell-centered data
461  if (!(this->_cell_centered_data.empty()))
462  {
463  std::map<std::string, const std::vector<Real> *>::iterator it = this->_cell_centered_data.begin();
464  const std::map<std::string, const std::vector<Real> *>::iterator end = this->_cell_centered_data.end();
465 
466  for (; it != end; ++it)
467  {
468  // write out the variable name, followed by a zero.
469  out_stream << (*it).first << " 0\n";
470 
471  const std::vector<Real> * the_array = (*it).second;
472 
473  // Loop over active elements, write out cell data. If second-order cells
474  // are split into sub-elements, the sub-elements inherit their parent's
475  // cell-centered data.
476  for (const auto & elem : mesh.active_element_ptr_range())
477  {
478  // Use the element's ID to find the value.
479  libmesh_assert_less (elem->id(), the_array->size());
480  const Real the_value = the_array->operator[](elem->id());
481 
482  if (this->subdivide_second_order())
483  for (unsigned int se=0; se < elem->n_sub_elem(); se++)
484  out_stream << the_value << " ";
485  else
486  out_stream << the_value << " ";
487  }
488 
489  out_stream << "\n\n";
490  }
491  }
492 
493 
494  // optionally write the data
495  if ((solution_names != libmesh_nullptr) && (v != libmesh_nullptr))
496  {
497  const unsigned int n_vars = solution_names->size();
498 
499  if (!(v->size() == mesh.n_nodes()*n_vars))
500  libMesh::err << "ERROR: v->size()=" << v->size()
501  << ", mesh.n_nodes()=" << mesh.n_nodes()
502  << ", n_vars=" << n_vars
503  << ", mesh.n_nodes()*n_vars=" << mesh.n_nodes()*n_vars
504  << "\n";
505 
506  libmesh_assert_equal_to (v->size(), mesh.n_nodes()*n_vars);
507 
508  for (unsigned int c=0; c<n_vars; c++)
509  {
510 
511 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
512 
513  // in case of complex data, write _three_ data sets
514  // for each component
515 
516  // this is the real part
517  out_stream << "r_" << (*solution_names)[c] << " 1\n";
518 
519  for (unsigned int n=0; n<mesh.n_nodes(); n++)
520  out_stream << (*v)[n*n_vars + c].real() << " ";
521 
522  out_stream << "\n\n";
523 
524  // this is the imaginary part
525  out_stream << "i_" << (*solution_names)[c] << " 1\n";
526 
527  for (unsigned int n=0; n<mesh.n_nodes(); n++)
528  out_stream << (*v)[n*n_vars + c].imag() << " ";
529 
530  out_stream << "\n\n";
531 
532  // this is the magnitude
533  out_stream << "a_" << (*solution_names)[c] << " 1\n";
534  for (unsigned int n=0; n<mesh.n_nodes(); n++)
535  out_stream << std::abs((*v)[n*n_vars + c]) << " ";
536 
537  out_stream << "\n\n";
538 
539 #else
540 
541  out_stream << (*solution_names)[c] << " 1\n";
542 
543  for (unsigned int n=0; n<mesh.n_nodes(); n++)
544  out_stream << (*v)[n*n_vars + c] << " ";
545 
546  out_stream << "\n\n";
547 
548 #endif
549  }
550 
551  }
552 
553  // If we wrote any variables, we have to close the variable section now
554  if (write_variable)
555  out_stream << "endvars\n";
556 
557 
558  // end of the file
559  out_stream << "\nendgmv\n";
560 
561 #endif
562 }
OStreamProxy err
double abs(double a)
virtual dof_id_type n_active_elem() const =0
virtual const Point & point(const dof_id_type i) const =0
unsigned int n_partitions() const
Definition: mesh_base.h:833
const class libmesh_nullptr_t libmesh_nullptr
unsigned int & ascii_precision()
Return/set the precision to use when writing ASCII files.
IterBase * end
Also have a polymorphic pointer to the end object, this prevents iterating past the end...
virtual SimpleRange< element_iterator > active_element_ptr_range()=0
const MT & mesh() const
Definition: mesh_output.h:216
This class defines an abstract interface for Mesh output.
Definition: mesh_output.h:53
const unsigned int n_vars
Definition: tecplot_io.C:68
bool & partitioning()
Flag indicating whether or not to write the partitioning information for the mesh.
Definition: gmv_io.h:107
long double max(long double a, double b)
This is the MeshBase class.
Definition: mesh_base.h:68
libmesh_assert(j)
void write_ascii_old_impl(const std::string &, const std::vector< Number > *=libmesh_nullptr, const std::vector< std::string > *=libmesh_nullptr)
This method implements writing a mesh with nodal data to a specified file where the nodal data and va...
Definition: gmv_io.C:569
bool & write_subdomain_id_as_material()
Flag to write element subdomain_id&#39;s as GMV "materials" instead of element processor_id&#39;s.
Definition: gmv_io.h:114
std::map< std::string, const std::vector< Real > * > _cell_centered_data
Storage for arbitrary cell-centered data.
Definition: gmv_io.h:225
bool & p_levels()
Flag indicating whether or not to write p level information for p refined meshes. ...
Definition: gmv_io.h:126
bool & subdivide_second_order()
Flag indicating whether or not to subdivide second order elements.
Definition: gmv_io.h:120
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
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 569 of file gmv_io.C.

References _cell_centered_data, std::abs(), libMesh::MeshBase::active_element_ptr_range(), libMesh::MeshOutput< MeshBase >::ascii_precision(), libMesh::Elem::build(), end, libMesh::Utility::enum_to_string(), libMesh::err, libMesh::FIRST, libMesh::Elem::first_order_equivalent_type(), libMesh::HEX20, libMesh::HEX27, libMesh::HEX8, 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(), n_vars, 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().

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

References _cell_centered_data, std::abs(), libMesh::MeshBase::active_element_ptr_range(), 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, p_levels(), partitioning(), libMesh::MeshBase::point(), libMesh::processor_id(), and write_subdomain_id_as_material().

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

1235 {
1236  // We use the file-scope global variable eletypes for mapping nodes
1237  // from GMV to libmesh indices, so initialize that data now.
1238  init_eletypes();
1239 
1240  // Get a reference to the mesh
1242 
1243  // This is a parallel_only function
1244  const dof_id_type n_active_elem = mesh.n_active_elem();
1245 
1247  return;
1248 
1249  std::ofstream out_stream (fname.c_str());
1250 
1251  libmesh_assert (out_stream.good());
1252 
1253  unsigned int mesh_max_p_level = 0;
1254 
1255  std::string buffer;
1256 
1257  // Begin interfacing with the GMV data file
1258  {
1259  // write the nodes
1260  buffer = "gmvinput";
1261  out_stream.write(buffer.c_str(), buffer.size());
1262 
1263  buffer = "ieeei4r4";
1264  out_stream.write(buffer.c_str(), buffer.size());
1265  }
1266 
1267 
1268 
1269  // write the nodes
1270  {
1271  buffer = "nodes ";
1272  out_stream.write(buffer.c_str(), buffer.size());
1273 
1274  unsigned int tempint = mesh.n_nodes();
1275  out_stream.write(reinterpret_cast<char *>(&tempint), sizeof(unsigned int));
1276 
1277  // write the x coordinate
1278  std::vector<float> temp(mesh.n_nodes());
1279  for (unsigned int v=0; v<mesh.n_nodes(); v++)
1280  temp[v] = static_cast<float>(mesh.point(v)(0));
1281  out_stream.write(reinterpret_cast<char *>(&temp[0]), sizeof(float)*mesh.n_nodes());
1282 
1283  // write the y coordinate
1284  for (unsigned int v=0; v<mesh.n_nodes(); v++)
1285  {
1286 #if LIBMESH_DIM > 1
1287  temp[v] = static_cast<float>(mesh.point(v)(1));
1288 #else
1289  temp[v] = 0.;
1290 #endif
1291  }
1292  out_stream.write(reinterpret_cast<char *>(&temp[0]), sizeof(float)*mesh.n_nodes());
1293 
1294  // write the z coordinate
1295  for (unsigned int v=0; v<mesh.n_nodes(); v++)
1296  {
1297 #if LIBMESH_DIM > 2
1298  temp[v] = static_cast<float>(mesh.point(v)(2));
1299 #else
1300  temp[v] = 0.;
1301 #endif
1302  }
1303  out_stream.write(reinterpret_cast<char *>(&temp[0]), sizeof(float)*mesh.n_nodes());
1304  }
1305 
1306 
1307  // write the connectivity
1308  {
1309  buffer = "cells ";
1310  out_stream.write(buffer.c_str(), buffer.size());
1311 
1312  unsigned int tempint = n_active_elem;
1313  out_stream.write(reinterpret_cast<char *>(&tempint), sizeof(unsigned int));
1314 
1315  for (const auto & elem : mesh.active_element_ptr_range())
1316  {
1317  mesh_max_p_level = std::max(mesh_max_p_level,
1318  elem->p_level());
1319 
1320  // The ElementDefinition label contains the GMV name
1321  // and the number of nodes for the respective
1322  // element.
1323  const ElementDefinition & ed = eletypes[elem->type()];
1324 
1325  // The ElementDefinition labels all have a name followed by a
1326  // whitespace and then the number of nodes. To write the
1327  // string in binary, we want to drop everything after that
1328  // space, and then pad the string out to 8 characters.
1329  buffer = ed.label;
1330 
1331  // Erase everything from the first whitespace character to the end of the string.
1332  buffer.erase(buffer.find_first_of(' '), std::string::npos);
1333 
1334  // Append whitespaces until the string is exactly 8 characters long.
1335  while (buffer.size() < 8)
1336  buffer.insert(buffer.end(), ' ');
1337 
1338  // Debugging:
1339  // libMesh::out << "Writing element with name = '" << buffer << "'" << std::endl;
1340 
1341  // Finally, write the 8 character stream to file.
1342  out_stream.write(buffer.c_str(), buffer.size());
1343 
1344  // Write the number of nodes that this type of element has.
1345  // Note: don't want to use the number of nodes of the element,
1346  // because certain elements, like QUAD9 and HEX27 only support
1347  // being written out as lower-order elements (QUAD8 and HEX20,
1348  // respectively).
1349  tempint = ed.node_map.size();
1350  out_stream.write(reinterpret_cast<char *>(&tempint), sizeof(unsigned int));
1351 
1352  // Write the element connectivity
1353  for (std::size_t i=0; i<ed.node_map.size(); i++)
1354  {
1355  dof_id_type id = elem->node_id(ed.node_map[i]) + 1;
1356  out_stream.write(reinterpret_cast<char *>(&id), sizeof(dof_id_type));
1357  }
1358  }
1359  }
1360 
1361 
1362 
1363  // optionally write the partition information
1364  if (this->partitioning())
1365  {
1366  if (this->write_subdomain_id_as_material())
1367  libmesh_error_msg("Not yet supported in GMVIO::write_binary");
1368 
1369  else
1370  {
1371  buffer = "material";
1372  out_stream.write(buffer.c_str(), buffer.size());
1373 
1374  unsigned int tmpint = mesh.n_processors();
1375  out_stream.write(reinterpret_cast<char *>(&tmpint), sizeof(unsigned int));
1376 
1377  tmpint = 0; // IDs are cell based
1378  out_stream.write(reinterpret_cast<char *>(&tmpint), sizeof(unsigned int));
1379 
1380  for (unsigned int proc=0; proc<mesh.n_processors(); proc++)
1381  {
1382  // We have to write exactly 8 bytes. This means that
1383  // the processor id can be up to 3 characters long, and
1384  // we pad it with blank characters on the end.
1385  std::ostringstream oss;
1386  oss << "proc_" << std::setw(3) << std::left << proc;
1387  out_stream.write(oss.str().c_str(), oss.str().size());
1388  }
1389 
1390  std::vector<unsigned int> proc_id (n_active_elem);
1391 
1392  unsigned int n=0;
1393 
1394  // We no longer write sub-elems in GMV, so just assign a
1395  // processor id value to each element.
1396  for (const auto & elem : mesh.active_element_ptr_range())
1397  proc_id[n++] = elem->processor_id() + 1;
1398 
1399  out_stream.write(reinterpret_cast<char *>(&proc_id[0]),
1400  sizeof(unsigned int)*proc_id.size());
1401  }
1402  }
1403 
1404  // If there are *any* variables at all in the system (including
1405  // p level, or arbitrary cell-based data)
1406  // to write, the gmv file needs to contain the word "variable"
1407  // on a line by itself.
1408  bool write_variable = false;
1409 
1410  // 1.) p-levels
1411  if (this->p_levels() && mesh_max_p_level)
1412  write_variable = true;
1413 
1414  // 2.) solution data
1415  if ((solution_names != libmesh_nullptr) && (vec != libmesh_nullptr))
1416  write_variable = true;
1417 
1418  // // 3.) cell-centered data - unsupported
1419  // if (!(this->_cell_centered_data.empty()))
1420  // write_variable = true;
1421 
1422  if (write_variable)
1423  {
1424  buffer = "variable";
1425  out_stream.write(buffer.c_str(), buffer.size());
1426  }
1427 
1428  // optionally write the partition information
1429  if (this->p_levels() && mesh_max_p_level)
1430  {
1431  unsigned int n_floats = n_active_elem;
1432  for (unsigned int i=0; i != mesh.mesh_dimension(); ++i)
1433  n_floats *= 2;
1434 
1435  std::vector<float> temp(n_floats);
1436 
1437  buffer = "p_level";
1438  out_stream.write(buffer.c_str(), buffer.size());
1439 
1440  unsigned int tempint = 0; // p levels are cell data
1441  out_stream.write(reinterpret_cast<char *>(&tempint), sizeof(unsigned int));
1442 
1443  unsigned int n=0;
1444  for (const auto & elem : mesh.active_element_ptr_range())
1445  for (unsigned int se=0; se<elem->n_sub_elem(); se++)
1446  temp[n++] = static_cast<float>( elem->p_level() );
1447 
1448  out_stream.write(reinterpret_cast<char *>(&temp[0]),
1449  sizeof(float)*n_floats);
1450  }
1451 
1452 
1453  // optionally write cell-centered data
1454  if (!(this->_cell_centered_data.empty()))
1455  libMesh::err << "Cell-centered data not (yet) supported in binary I/O mode!" << std::endl;
1456 
1457 
1458 
1459 
1460  // optionally write the data
1461  if ((solution_names != libmesh_nullptr) &&
1462  (vec != libmesh_nullptr))
1463  {
1464  std::vector<float> temp(mesh.n_nodes());
1465 
1466  const unsigned int n_vars =
1467  cast_int<unsigned int>(solution_names->size());
1468 
1469  for (unsigned int c=0; c<n_vars; c++)
1470  {
1471 
1472 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
1473  // for complex data, write three datasets
1474 
1475 
1476  // Real part
1477  buffer = "r_";
1478  out_stream.write(buffer.c_str(), buffer.size());
1479 
1480  buffer = (*solution_names)[c];
1481  out_stream.write(buffer.c_str(), buffer.size());
1482 
1483  unsigned int tempint = 1; // always do nodal data
1484  out_stream.write(reinterpret_cast<char *>(&tempint), sizeof(unsigned int));
1485 
1486  for (unsigned int n=0; n<mesh.n_nodes(); n++)
1487  temp[n] = static_cast<float>( (*vec)[n*n_vars + c].real() );
1488 
1489  out_stream.write(reinterpret_cast<char *>(&temp[0]), sizeof(float)*mesh.n_nodes());
1490 
1491 
1492  // imaginary part
1493  buffer = "i_";
1494  out_stream.write(buffer.c_str(), buffer.size());
1495 
1496  buffer = (*solution_names)[c];
1497  out_stream.write(buffer.c_str(), buffer.size());
1498 
1499  out_stream.write(reinterpret_cast<char *>(&tempint), sizeof(unsigned int));
1500 
1501  for (unsigned int n=0; n<mesh.n_nodes(); n++)
1502  temp[n] = static_cast<float>( (*vec)[n*n_vars + c].imag() );
1503 
1504  out_stream.write(reinterpret_cast<char *>(&temp[0]), sizeof(float)*mesh.n_nodes());
1505 
1506  // magnitude
1507  buffer = "a_";
1508  out_stream.write(buffer.c_str(), buffer.size());
1509  buffer = (*solution_names)[c];
1510  out_stream.write(buffer.c_str(), buffer.size());
1511 
1512  out_stream.write(reinterpret_cast<char *>(&tempint), sizeof(unsigned int));
1513 
1514  for (unsigned int n=0; n<mesh.n_nodes(); n++)
1515  temp[n] = static_cast<float>(std::abs((*vec)[n*n_vars + c]));
1516 
1517  out_stream.write(reinterpret_cast<char *>(&temp[0]), sizeof(float)*mesh.n_nodes());
1518 
1519 #else
1520 
1521  buffer = (*solution_names)[c];
1522  out_stream.write(buffer.c_str(), buffer.size());
1523 
1524  unsigned int tempint = 1; // always do nodal data
1525  out_stream.write(reinterpret_cast<char *>(&tempint), sizeof(unsigned int));
1526 
1527  for (unsigned int n=0; n<mesh.n_nodes(); n++)
1528  temp[n] = static_cast<float>((*vec)[n*n_vars + c]);
1529 
1530  out_stream.write(reinterpret_cast<char *>(&temp[0]), sizeof(float)*mesh.n_nodes());
1531 
1532 #endif
1533  }
1534  }
1535 
1536  // If we wrote any variables, we have to close the variable section now
1537  if (write_variable)
1538  {
1539  buffer = "endvars ";
1540  out_stream.write(buffer.c_str(), buffer.size());
1541  }
1542 
1543  // end the file
1544  buffer = "endgmv ";
1545  out_stream.write(buffer.c_str(), buffer.size());
1546 }
OStreamProxy err
double abs(double a)
virtual dof_id_type n_active_elem() const =0
virtual const Point & point(const dof_id_type i) const =0
processor_id_type n_processors() const
const class libmesh_nullptr_t libmesh_nullptr
virtual SimpleRange< element_iterator > active_element_ptr_range()=0
const MT & mesh() const
Definition: mesh_output.h:216
This class defines an abstract interface for Mesh output.
Definition: mesh_output.h:53
const unsigned int n_vars
Definition: tecplot_io.C:68
bool & partitioning()
Flag indicating whether or not to write the partitioning information for the mesh.
Definition: gmv_io.h:107
long double max(long double a, double b)
This is the MeshBase class.
Definition: mesh_base.h:68
libmesh_assert(j)
bool & write_subdomain_id_as_material()
Flag to write element subdomain_id&#39;s as GMV "materials" instead of element processor_id&#39;s.
Definition: gmv_io.h:114
std::map< std::string, const std::vector< Real > * > _cell_centered_data
Storage for arbitrary cell-centered data.
Definition: gmv_io.h:225
bool & p_levels()
Flag indicating whether or not to write p level information for p refined meshes. ...
Definition: gmv_io.h:126
unsigned int mesh_dimension() const
Definition: mesh_base.C:148
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
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 1556 of file gmv_io.C.

References _cell_centered_data, _write_subdomain_id_as_material, std::abs(), libMesh::MeshBase::active_element_ptr_range(), libMesh::EquationSystems::build_discontinuous_solution_vector(), libMesh::EquationSystems::build_variable_names(), libMesh::EDGE2, libMesh::EDGE3, libMesh::EDGE4, 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().

1560 {
1561  std::vector<std::string> solution_names;
1562  std::vector<Number> v;
1563 
1564  // Get a reference to the mesh
1566 
1567  es.build_variable_names (solution_names, libmesh_nullptr, system_names);
1568  es.build_discontinuous_solution_vector (v, system_names);
1569 
1570  // These are parallel_only functions
1571  const unsigned int n_active_elem = mesh.n_active_elem();
1572 
1573  if (mesh.processor_id() != 0)
1574  return;
1575 
1576  std::ofstream out_stream(name.c_str());
1577 
1578  libmesh_assert (out_stream.good());
1579 
1580  // Begin interfacing with the GMV data file
1581  {
1582 
1583  // write the nodes
1584  out_stream << "gmvinput ascii" << std::endl << std::endl;
1585 
1586  // Compute the total weight
1587  {
1588  unsigned int tw=0;
1589 
1590  for (const auto & elem : mesh.active_element_ptr_range())
1591  tw += elem->n_nodes();
1592 
1593  out_stream << "nodes " << tw << std::endl;
1594  }
1595 
1596 
1597 
1598  // Write all the x values
1599  {
1600  for (const auto & elem : mesh.active_element_ptr_range())
1601  for (unsigned int n=0; n<elem->n_nodes(); n++)
1602  out_stream << elem->point(n)(0) << " ";
1603 
1604  out_stream << std::endl;
1605  }
1606 
1607 
1608  // Write all the y values
1609  {
1610  for (const auto & elem : mesh.active_element_ptr_range())
1611  for (unsigned int n=0; n<elem->n_nodes(); n++)
1612 #if LIBMESH_DIM > 1
1613  out_stream << elem->point(n)(1) << " ";
1614 #else
1615  out_stream << 0. << " ";
1616 #endif
1617 
1618  out_stream << std::endl;
1619  }
1620 
1621 
1622  // Write all the z values
1623  {
1624  for (const auto & elem : mesh.active_element_ptr_range())
1625  for (unsigned int n=0; n<elem->n_nodes(); n++)
1626 #if LIBMESH_DIM > 2
1627  out_stream << elem->point(n)(2) << " ";
1628 #else
1629  out_stream << 0. << " ";
1630 #endif
1631 
1632  out_stream << std::endl << std::endl;
1633  }
1634  }
1635 
1636 
1637 
1638  {
1639  // write the connectivity
1640 
1641  out_stream << "cells " << n_active_elem << std::endl;
1642 
1643  unsigned int nn=1;
1644 
1645  switch (mesh.mesh_dimension())
1646  {
1647  case 1:
1648  {
1649  for (const auto & elem : mesh.active_element_ptr_range())
1650  for (unsigned int se=0; se<elem->n_sub_elem(); se++)
1651  {
1652  if ((elem->type() == EDGE2) ||
1653  (elem->type() == EDGE3) ||
1654  (elem->type() == EDGE4)
1655 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1656  || (elem->type() == INFEDGE2)
1657 #endif
1658  )
1659  {
1660  out_stream << "line 2" << std::endl;
1661  for (unsigned int i=0; i<elem->n_nodes(); i++)
1662  out_stream << nn++ << " ";
1663 
1664  }
1665  else
1666  libmesh_error_msg("Unsupported 1D element type: " << Utility::enum_to_string(elem->type()));
1667 
1668  out_stream << std::endl;
1669  }
1670 
1671  break;
1672  }
1673 
1674  case 2:
1675  {
1676  for (const auto & elem : mesh.active_element_ptr_range())
1677  for (unsigned int se=0; se<elem->n_sub_elem(); se++)
1678  {
1679  if ((elem->type() == QUAD4) ||
1680  (elem->type() == QUAD8) || // Note: QUAD8 will be output as one central quad and
1681  // four surrounding triangles (though they will be written
1682  // to GMV as QUAD4s).
1683  (elem->type() == QUAD9)
1684 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1685  || (elem->type() == INFQUAD4)
1686  || (elem->type() == INFQUAD6)
1687 #endif
1688  )
1689  {
1690  out_stream << "quad 4" << std::endl;
1691  for (unsigned int i=0; i<elem->n_nodes(); i++)
1692  out_stream << nn++ << " ";
1693 
1694  }
1695  else if ((elem->type() == TRI3) ||
1696  (elem->type() == TRI6))
1697  {
1698  out_stream << "tri 3" << std::endl;
1699  for (unsigned int i=0; i<elem->n_nodes(); i++)
1700  out_stream << nn++ << " ";
1701 
1702  }
1703  else
1704  libmesh_error_msg("Unsupported 2D element type: " << Utility::enum_to_string(elem->type()));
1705 
1706  out_stream << std::endl;
1707  }
1708 
1709  break;
1710  }
1711 
1712 
1713  case 3:
1714  {
1715  for (const auto & elem : mesh.active_element_ptr_range())
1716  for (unsigned int se=0; se<elem->n_sub_elem(); se++)
1717  {
1718  if ((elem->type() == HEX8) ||
1719  (elem->type() == HEX20) ||
1720  (elem->type() == HEX27)
1721 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1722  || (elem->type() == INFHEX8)
1723  || (elem->type() == INFHEX16)
1724  || (elem->type() == INFHEX18)
1725 #endif
1726  )
1727  {
1728  out_stream << "phex8 8" << std::endl;
1729  for (unsigned int i=0; i<elem->n_nodes(); i++)
1730  out_stream << nn++ << " ";
1731  }
1732  else if ((elem->type() == PRISM6) ||
1733  (elem->type() == PRISM15) ||
1734  (elem->type() == PRISM18)
1735 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1736  || (elem->type() == INFPRISM6)
1737  || (elem->type() == INFPRISM12)
1738 #endif
1739  )
1740  {
1741  out_stream << "pprism6 6" << std::endl;
1742  for (unsigned int i=0; i<elem->n_nodes(); i++)
1743  out_stream << nn++ << " ";
1744  }
1745  else if ((elem->type() == TET4) ||
1746  (elem->type() == TET10))
1747  {
1748  out_stream << "tet 4" << std::endl;
1749  for (unsigned int i=0; i<elem->n_nodes(); i++)
1750  out_stream << nn++ << " ";
1751  }
1752  else
1753  libmesh_error_msg("Unsupported 3D element type: " << Utility::enum_to_string(elem->type()));
1754 
1755  out_stream << std::endl;
1756  }
1757 
1758  break;
1759  }
1760 
1761  default:
1762  libmesh_error_msg("Unsupported mesh dimension: " << mesh.mesh_dimension());
1763  }
1764 
1765  out_stream << std::endl;
1766  }
1767 
1768 
1769 
1770  // optionally write the partition information
1771  if (write_partitioning)
1772  {
1774  libmesh_error_msg("Not yet supported in GMVIO::write_discontinuous_gmv");
1775 
1776  else
1777  {
1778  out_stream << "material "
1779  << mesh.n_processors()
1780  << " 0"<< std::endl;
1781 
1782  for (unsigned int proc=0; proc<mesh.n_processors(); proc++)
1783  out_stream << "proc_" << proc << std::endl;
1784 
1785  for (const auto & elem : mesh.active_element_ptr_range())
1786  out_stream << elem->processor_id()+1 << std::endl;
1787 
1788  out_stream << std::endl;
1789  }
1790  }
1791 
1792 
1793  // Writing cell-centered data is not yet supported in discontinuous GMV files.
1794  if (!(this->_cell_centered_data.empty()))
1795  {
1796  libMesh::err << "Cell-centered data not (yet) supported for discontinuous GMV files!" << std::endl;
1797  }
1798 
1799 
1800 
1801  // write the data
1802  {
1803  const unsigned int n_vars =
1804  cast_int<unsigned int>(solution_names.size());
1805 
1806  // libmesh_assert_equal_to (v.size(), tw*n_vars);
1807 
1808  out_stream << "variable" << std::endl;
1809 
1810 
1811  for (unsigned int c=0; c<n_vars; c++)
1812  {
1813 
1814 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
1815 
1816  // in case of complex data, write _tree_ data sets
1817  // for each component
1818 
1819  // this is the real part
1820  out_stream << "r_" << solution_names[c] << " 1" << std::endl;
1821  for (const auto & elem : mesh.active_element_ptr_range())
1822  for (unsigned int n=0; n<elem->n_nodes(); n++)
1823  out_stream << v[(n++)*n_vars + c].real() << " ";
1824  out_stream << std::endl << std::endl;
1825 
1826 
1827  // this is the imaginary part
1828  out_stream << "i_" << solution_names[c] << " 1" << std::endl;
1829  for (const auto & elem : mesh.active_element_ptr_range())
1830  for (unsigned int n=0; n<elem->n_nodes(); n++)
1831  out_stream << v[(n++)*n_vars + c].imag() << " ";
1832  out_stream << std::endl << std::endl;
1833 
1834  // this is the magnitude
1835  out_stream << "a_" << solution_names[c] << " 1" << std::endl;
1836  for (const auto & elem : mesh.active_element_ptr_range())
1837  for (unsigned int n=0; n<elem->n_nodes(); n++)
1838  out_stream << std::abs(v[(n++)*n_vars + c]) << " ";
1839  out_stream << std::endl << std::endl;
1840 
1841 #else
1842 
1843  out_stream << solution_names[c] << " 1" << std::endl;
1844  {
1845  unsigned int nn=0;
1846 
1847  for (const auto & elem : mesh.active_element_ptr_range())
1848  for (unsigned int n=0; n<elem->n_nodes(); n++)
1849  out_stream << v[(nn++)*n_vars + c] << " ";
1850  }
1851  out_stream << std::endl << std::endl;
1852 
1853 #endif
1854 
1855  }
1856 
1857  out_stream << "endvars" << std::endl;
1858  }
1859 
1860 
1861  // end of the file
1862  out_stream << std::endl << "endgmv" << std::endl;
1863 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:39
OStreamProxy err
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
virtual SimpleRange< element_iterator > active_element_ptr_range()=0
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
Fill the input vector soln with solution values.
This is the MeshBase class.
Definition: mesh_base.h:68
libmesh_assert(j)
std::map< std::string, const std::vector< Real > * > _cell_centered_data
Storage for arbitrary cell-centered data.
Definition: gmv_io.h:225
bool _write_subdomain_id_as_material
Flag to write element subdomain_id&#39;s as GMV "materials" instead of element processor_id&#39;s.
Definition: gmv_io.h:207
std::string enum_to_string(const T e)
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
Fill the input vector var_names with the names of the variables for each system.
unsigned int mesh_dimension() const
Definition: mesh_base.C:148
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()
Flag indicating whether or not to write a binary file.
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)
This method implements writing a mesh with nodal data to a specified file where the nodal data and va...
Definition: gmv_io.C:569
void write_binary(const std::string &, const std::vector< Number > *=libmesh_nullptr, const std::vector< std::string > *=libmesh_nullptr)
This method implements writing a mesh with nodal data to a specified file where the nodal data and va...
Definition: gmv_io.C:1232
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 ( )

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
Flag to write element subdomain_id&#39;s as GMV "materials" instead of element processor_id&#39;s.
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().

std::vector<bool> libMesh::MeshInput< MeshBase >::elems_of_dimension
protectedinherited

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