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

This class implements writing meshes and solutions in Ensight's Gold format. More...

#include <ensight_io.h>

Inheritance diagram for libMesh::EnsightIO:
[legend]

Classes

struct  Scalars
 
struct  SystemVars
 
struct  Vectors
 

Public Member Functions

 EnsightIO (const std::string &filename, const EquationSystems &eq)
 Constructor. More...
 
 ~EnsightIO ()
 Empty destructor. More...
 
void add_scalar (const std::string &system, const std::string &scalar_description, const std::string &s)
 Tell the EnsightIO interface to output the finite element (not SCALAR) variable named "s". More...
 
void add_vector (const std::string &system, const std::string &vec_description, const std::string &u, const std::string &v)
 Tell the EnsightIO interface that the variables (u,v) constitute a vector. More...
 
void add_vector (const std::string &system, const std::string &vec_description, const std::string &u, const std::string &v, const std::string &w)
 Tell the EnsightIO interface that the variables (u, v, w) constitute a vector. More...
 
void write (Real time=0)
 Calls write_ascii() and write_case(). More...
 
virtual void write (const std::string &name) libmesh_override
 Calls this->write(0);. 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 std::vector< Number > &, const std::vector< std::string > &)
 This method implements writing a mesh with nodal data to a specified file where the nodal data and variable names are provided. 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

const MeshBasemesh () const
 

Protected Attributes

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 Types

typedef std::map< std::string, SystemVarssystem_vars_map_t
 

Private Member Functions

void write_ascii (Real time=0)
 
void write_scalar_ascii (const std::string &sys, const std::string &var)
 
void write_vector_ascii (const std::string &sys, const std::vector< std::string > &vec, const std::string &var_name)
 
void write_solution_ascii ()
 
void write_geometry_ascii ()
 
void write_case ()
 

Static Private Member Functions

static std::map< ElemType, std::string > build_element_map ()
 

Private Attributes

std::string _ensight_file_name
 
std::vector< Real_time_steps
 
system_vars_map_t _system_vars_map
 
const EquationSystems_equation_systems
 

Static Private Attributes

static std::map< ElemType, std::string > _element_map = EnsightIO::build_element_map()
 

Detailed Description

This class implements writing meshes and solutions in Ensight's Gold format.

Author
Camata
Date
2009
Author
J. W. Peterson (refactoring and iostreams implementation)
Date
2016

Definition at line 47 of file ensight_io.h.

Member Typedef Documentation

typedef std::map<std::string, SystemVars> libMesh::EnsightIO::system_vars_map_t
private

Definition at line 148 of file ensight_io.h.

Constructor & Destructor Documentation

libMesh::EnsightIO::EnsightIO ( const std::string &  filename,
const EquationSystems eq 
)

Constructor.

Definition at line 60 of file ensight_io.C.

References _ensight_file_name, _equation_systems, libMesh::ParallelObject::n_processors(), and libMesh::ParallelObject::processor_id().

61  :
62  MeshOutput<MeshBase> (eq.get_mesh()),
64 {
66  _ensight_file_name = filename;
67  else
68  {
69  std::ostringstream tmp_file;
70  tmp_file << filename << "_rank" << _equation_systems.processor_id();
71  _ensight_file_name = tmp_file.str();
72  }
73 }
std::string _ensight_file_name
Definition: ensight_io.h:144
processor_id_type n_processors() const
PetscBool eq
const EquationSystems & _equation_systems
Definition: ensight_io.h:152
processor_id_type processor_id() const
libMesh::EnsightIO::~EnsightIO ( )

Empty destructor.

Definition at line 60 of file ensight_io.h.

References add_scalar(), add_vector(), libMesh::Quality::name(), libMesh::Real, and write().

60 {}

Member Function Documentation

void libMesh::EnsightIO::add_scalar ( const std::string &  system,
const std::string &  scalar_description,
const std::string &  s 
)

Tell the EnsightIO interface to output the finite element (not SCALAR) variable named "s".

Note
You must call add_scalar() or add_vector() (see below) at least once, otherwise only the Mesh will be written out.

Definition at line 117 of file ensight_io.C.

References _equation_systems, _system_vars_map, libMesh::EnsightIO::Scalars::description, libMesh::EquationSystems::get_system(), libMesh::EquationSystems::has_system(), and libMesh::libmesh_assert().

Referenced by ~EnsightIO().

120 {
122  libmesh_assert(_equation_systems.get_system(system_name).has_variable(s));
123 
124  Scalars scl;
125  scl.description = scl_description;
126  scl.scalar_name = s;
127 
128  _system_vars_map[system_name].EnsightScalars.push_back(scl);
129 }
bool has_system(const std::string &name) const
libmesh_assert(j)
system_vars_map_t _system_vars_map
Definition: ensight_io.h:149
const EquationSystems & _equation_systems
Definition: ensight_io.h:152
const T_sys & get_system(const std::string &name) const
void libMesh::EnsightIO::add_vector ( const std::string &  system,
const std::string &  vec_description,
const std::string &  u,
const std::string &  v 
)

Tell the EnsightIO interface that the variables (u,v) constitute a vector.

Note
u and v must have the same FEType, and be defined in the same system.

Definition at line 77 of file ensight_io.C.

References _equation_systems, _system_vars_map, libMesh::EnsightIO::Vectors::description, libMesh::EquationSystems::get_system(), libMesh::EquationSystems::has_system(), and libMesh::libmesh_assert().

Referenced by ~EnsightIO().

81 {
83  libmesh_assert (_equation_systems.get_system(system_name).has_variable(u));
84  libmesh_assert (_equation_systems.get_system(system_name).has_variable(v));
85 
86  Vectors vec;
87  vec.description = vec_description;
88  vec.components.push_back(u);
89  vec.components.push_back(v);
90 
91  _system_vars_map[system_name].EnsightVectors.push_back(vec);
92 }
bool has_system(const std::string &name) const
libmesh_assert(j)
system_vars_map_t _system_vars_map
Definition: ensight_io.h:149
const EquationSystems & _equation_systems
Definition: ensight_io.h:152
const T_sys & get_system(const std::string &name) const
void libMesh::EnsightIO::add_vector ( const std::string &  system,
const std::string &  vec_description,
const std::string &  u,
const std::string &  v,
const std::string &  w 
)

Tell the EnsightIO interface that the variables (u, v, w) constitute a vector.

Note
Requires a 3D mesh, u, v, and w must have the same FEType, and must be defined in the same system.

Definition at line 96 of file ensight_io.C.

References _equation_systems, _system_vars_map, libMesh::EnsightIO::Vectors::description, libMesh::EquationSystems::get_system(), libMesh::EquationSystems::has_system(), and libMesh::libmesh_assert().

101 {
103  libmesh_assert(_equation_systems.get_system(system_name).has_variable(u));
104  libmesh_assert(_equation_systems.get_system(system_name).has_variable(v));
105  libmesh_assert(_equation_systems.get_system(system_name).has_variable(w));
106 
107  Vectors vec;
108  vec.description = vec_name;
109  vec.components.push_back(u);
110  vec.components.push_back(v);
111  vec.components.push_back(w);
112  _system_vars_map[system_name].EnsightVectors.push_back(vec);
113 }
bool has_system(const std::string &name) const
libmesh_assert(j)
system_vars_map_t _system_vars_map
Definition: ensight_io.h:149
const EquationSystems & _equation_systems
Definition: ensight_io.h:152
const T_sys & get_system(const std::string &name) const
unsigned int& libMesh::MeshOutput< MeshBase >::ascii_precision ( )
inherited

Return/set the precision to use when writing ASCII files.

By default we use numeric_limits<Real>::digits10 + 2, which should be enough to write out to ASCII and get the exact same Real back when reading in.

Referenced by libMesh::TecplotIO::write_ascii(), libMesh::GMVIO::write_ascii_new_impl(), and libMesh::GMVIO::write_ascii_old_impl().

std::map< ElemType, std::string > libMesh::EnsightIO::build_element_map ( )
staticprivate

Definition at line 40 of file ensight_io.C.

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

41 {
42  std::map<ElemType, std::string> ret;
43  ret[EDGE2] = "bar2";
44  ret[EDGE3] = "bar3";
45  ret[QUAD4] = "quad4";
46  ret[QUAD8] = "quad8";
47  // ret[QUAD9] = "quad9"; // not supported
48  ret[TRI3] = "tria3";
49  ret[TRI6] = "tria6";
50  ret[TET4] = "tetra4";
51  ret[TET10] = "tetra10";
52  ret[HEX8] = "hexa8";
53  ret[HEX20] = "hexa20";
54  // ret[HEX27] = "HEX27"; // not supported
55  ret[PYRAMID5] = "pyramid5";
56  return ret;
57 }
const MeshBase & libMesh::MeshOutput< MeshBase >::mesh ( ) const
protectedinherited
void libMesh::EnsightIO::write ( Real  time = 0)

Calls write_ascii() and write_case().

Writes case, mesh, and solution files named: name.case (contains a description of other files) name.geo000 (mesh) name_{varname}.scl000 (one file per scalar variable) name_{vecname}.vec000 (one file per vector variable)

Definition at line 147 of file ensight_io.C.

References write_ascii(), and write_case().

Referenced by write(), and ~EnsightIO().

148 {
149  this->write_ascii(time);
150  this->write_case();
151 }
void write_ascii(Real time=0)
Definition: ensight_io.C:155
void libMesh::EnsightIO::write ( const std::string &  name)
virtual

Calls this->write(0);.

Implements libMesh::MeshOutput< MeshBase >.

Definition at line 135 of file ensight_io.C.

References _ensight_file_name, libMesh::MeshOutput< MeshBase >::_is_parallel_format, libMesh::MeshOutput< MeshBase >::mesh(), libMesh::Quality::name(), and write().

136 {
137  // We may need to gather a DistributedMesh to output it, making that
138  // const qualifier in our constructor a dirty lie
139  MeshSerializer serialize(const_cast<MeshBase &>(this->mesh()), !_is_parallel_format);
140 
142  this->write();
143 }
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:39
std::string _ensight_file_name
Definition: ensight_io.h:144
const bool _is_parallel_format
Flag specifying whether this format is parallel-capable.
Definition: mesh_output.h:141
const MeshBase & mesh() const
void write(Real time=0)
Calls write_ascii() and write_case().
Definition: ensight_io.C:147
void libMesh::EnsightIO::write_ascii ( Real  time = 0)
private

Definition at line 155 of file ensight_io.C.

References _time_steps, write_geometry_ascii(), and write_solution_ascii().

Referenced by write().

156 {
157  _time_steps.push_back(time);
158 
159  this->write_geometry_ascii();
160  this->write_solution_ascii();
161 }
void write_solution_ascii()
Definition: ensight_io.C:342
void write_geometry_ascii()
Definition: ensight_io.C:165
std::vector< Real > _time_steps
Definition: ensight_io.h:145
void libMesh::EnsightIO::write_case ( )
private

Definition at line 287 of file ensight_io.C.

References _ensight_file_name, _system_vars_map, _time_steps, libMesh::EnsightIO::Vectors::description, libMesh::EnsightIO::Scalars::description, and libMesh::EnsightIO::Scalars::scalar_name.

Referenced by write().

288 {
289  std::ostringstream case_file;
290  case_file << _ensight_file_name << ".case";
291 
292  // Open a stream for writing the case file.
293  std::ofstream case_stream(case_file.str().c_str());
294 
295  case_stream << "FORMAT\n";
296  case_stream << "type: ensight gold\n\n";
297  case_stream << "GEOMETRY\n";
298  case_stream << "model: 1 " << _ensight_file_name << ".geo" << "***\n";
299 
300  system_vars_map_t::iterator sys_it = _system_vars_map.begin();
301  const system_vars_map_t::iterator sys_end = _system_vars_map.end();
302 
303  // Write Variable per node section
304  if (sys_it != sys_end)
305  case_stream << "\n\nVARIABLE\n";
306 
307  for (; sys_it != sys_end; ++sys_it)
308  {
309  for (std::size_t i=0; i < sys_it->second.EnsightScalars.size(); i++)
310  {
311  Scalars scalar = sys_it->second.EnsightScalars[i];
312  case_stream << "scalar per node: 1 "
313  << scalar.description << " "
314  << _ensight_file_name << "_" << scalar.scalar_name << ".scl***\n";
315  }
316 
317  for (std::size_t i=0; i < sys_it->second.EnsightVectors.size(); i++)
318  {
319  Vectors vec = sys_it->second.EnsightVectors[i];
320  case_stream << "vector per node: 1 "
321  << vec.description << " "
322  << _ensight_file_name << "_" << vec.description << ".vec***\n";
323  }
324 
325  // Write time step section
326  if (_time_steps.size() != 0)
327  {
328  case_stream << "\n\nTIME\n";
329  case_stream << "time set: 1\n";
330  case_stream << "number of steps: " << std::setw(10) << _time_steps.size() << "\n";
331  case_stream << "filename start number: " << std::setw(10) << 0 << "\n";
332  case_stream << "filename increment: " << std::setw(10) << 1 << "\n";
333  case_stream << "time values:\n";
334  for (std::size_t i = 0; i < _time_steps.size(); i++)
335  case_stream << std::setw(12) << std::setprecision(5) << std::scientific << _time_steps[i] << "\n";
336  }
337  }
338 }
std::string _ensight_file_name
Definition: ensight_io.h:144
system_vars_map_t _system_vars_map
Definition: ensight_io.h:149
std::vector< Real > _time_steps
Definition: ensight_io.h:145
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::EnsightIO::write_geometry_ascii ( )
private

Definition at line 165 of file ensight_io.C.

References _element_map, _ensight_file_name, _time_steps, libMesh::MeshBase::active_local_element_ptr_range(), libMesh::HEX27, libMesh::MeshOutput< MT >::mesh(), and libMesh::QUAD9.

Referenced by write_ascii().

166 {
167  std::ostringstream file;
168  file << _ensight_file_name
169  << ".geo"
170  << std::setw(3)
171  << std::setprecision(0)
172  << std::setfill('0')
173  << std::right
174  << _time_steps.size()-1;
175 
176  // Open a stream to write the mesh
177  std::ofstream mesh_stream(file.str().c_str());
178 
179  mesh_stream << "EnSight Gold Geometry File Format\n";
180  mesh_stream << "Generated by \n";
181  mesh_stream << "node id off\n";
182  mesh_stream << "element id given\n";
183  mesh_stream << "part\n";
184  mesh_stream << std::setw(10) << 1 << "\n";
185  mesh_stream << "uns-elements\n";
186  mesh_stream << "coordinates\n";
187 
188  // mapping between nodal index and your coordinates
189  typedef std::map<int, Point> mesh_nodes_map_t;
190  typedef mesh_nodes_map_t::iterator mesh_nodes_iterator;
191  mesh_nodes_map_t mesh_nodes_map;
192 
193  // Map for grouping elements of the same type
194  typedef std::map<ElemType, std::vector<const Elem *>> ensight_parts_map_t;
195  typedef ensight_parts_map_t::iterator ensight_parts_iterator;
196  ensight_parts_map_t ensight_parts_map;
197 
198  const MeshBase & the_mesh = MeshOutput<MeshBase>::mesh();
199 
200  // Construct the various required maps
201  for (const auto & elem : the_mesh.active_local_element_ptr_range())
202  {
203  ensight_parts_map[elem->type()].push_back(elem);
204 
205  for (unsigned int i = 0; i < elem->n_nodes(); i++)
206  mesh_nodes_map[elem->node_id(i)] = elem->point(i);
207  }
208 
209  // Write number of local points
210  mesh_stream << std::setw(10) << mesh_nodes_map.size() << "\n";
211 
212  // write x, y, and z node positions, build mapping between
213  // ensight and libmesh node numbers.
214  std::map <int, int> ensight_node_index;
215  {
216  mesh_nodes_iterator no_it = mesh_nodes_map.begin();
217  const mesh_nodes_iterator no_end_it = mesh_nodes_map.end();
218 
219  for (unsigned direction=0; direction<3; ++direction)
220  {
221  for (int i = 1; no_it != no_end_it; ++no_it, i++)
222  {
223  mesh_stream << std::setw(12)
224  << std::setprecision(5)
225  << std::scientific
226  << no_it->second(direction)
227  << "\n";
228  ensight_node_index[no_it->first] = i;
229  }
230 
231  // Reset iterator to the beginning of the map
232  no_it = mesh_nodes_map.begin();
233  }
234  }
235 
236  // Write parts
237  {
238  ensight_parts_iterator parts_it = ensight_parts_map.begin();
239  const ensight_parts_iterator end_parts_it = ensight_parts_map.end();
240 
241  for (; parts_it != end_parts_it; ++parts_it)
242  {
243  // Look up this ElemType in the map, error if not present.
244  std::map<ElemType, std::string>::iterator name_it = _element_map.find(parts_it->first);
245  if (name_it == _element_map.end())
246  libmesh_error_msg("Error: Unsupported ElemType " << parts_it->first << " for EnsightIO.");
247 
248  // Write element type
249  mesh_stream << "\n" << name_it->second << "\n";
250 
251  std::vector<const Elem *> elem_ref = parts_it->second;
252 
253  // Write number of element
254  mesh_stream << std::setw(10) << elem_ref.size() << "\n";
255 
256  // Write element id
257  for (std::size_t i = 0; i < elem_ref.size(); i++)
258  mesh_stream << std::setw(10) << elem_ref[i]->id() << "\n";
259 
260  // Write connectivity
261  for (std::size_t i = 0; i < elem_ref.size(); i++)
262  {
263  for (unsigned int j = 0; j < elem_ref[i]->n_nodes(); j++)
264  {
265  // tests!
266  if (parts_it->first == QUAD9 && i==4)
267  continue;
268 
269  // tests!
270  if (parts_it->first == HEX27 &&
271  (i==4 || i ==10 || i == 12 ||
272  i == 13 || i ==14 || i == 16 || i == 22))
273  continue;
274 
275  mesh_stream << std::setw(10) << ensight_node_index[elem_ref[i]->node_id(j)];
276  }
277  mesh_stream << "\n";
278  }
279  }
280  }
281 }
std::string _ensight_file_name
Definition: ensight_io.h:144
const MT & mesh() const
Definition: mesh_output.h:216
static std::map< ElemType, std::string > _element_map
Definition: ensight_io.h:155
std::vector< Real > _time_steps
Definition: ensight_io.h:145
virtual void libMesh::MeshOutput< MeshBase >::write_nodal_data ( const std::string &  ,
const std::vector< Number > &  ,
const std::vector< std::string > &   
)
virtualinherited

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

Reimplemented in libMesh::ExodusII_IO, libMesh::NameBasedIO, libMesh::GmshIO, libMesh::Nemesis_IO, libMesh::VTKIO, libMesh::UCDIO, libMesh::GMVIO, libMesh::MEDITIO, libMesh::GnuPlotIO, and libMesh::TecplotIO.

Definition at line 96 of file mesh_output.h.

References libMesh::MeshOutput< MT >::ascii_precision(), libMesh::MeshOutput< MT >::mesh(), and libMesh::MeshOutput< MT >::write_nodal_data().

99  { libmesh_not_implemented(); }
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.

void libMesh::EnsightIO::write_scalar_ascii ( const std::string &  sys,
const std::string &  var 
)
private

Definition at line 361 of file ensight_io.C.

References _ensight_file_name, _equation_systems, _time_steps, libMesh::MeshBase::active_local_element_ptr_range(), libMesh::System::current_solution(), dim, libMesh::DofMap::dof_indices(), libMesh::System::get_dof_map(), libMesh::EquationSystems::get_system(), libMesh::libmesh_real(), libMesh::MeshOutput< MT >::mesh(), libMesh::MeshBase::mesh_dimension(), libMesh::FEInterface::nodal_soln(), libMesh::System::variable_number(), and libMesh::System::variable_type().

Referenced by write_solution_ascii().

363 {
364  // Construct scalar variable filename
365  std::ostringstream scl_file;
366  scl_file << _ensight_file_name
367  << "_"
368  << var_name
369  << ".scl"
370  << std::setw(3)
371  << std::setprecision(0)
372  << std::setfill('0')
373  << std::right
374  << _time_steps.size()-1;
375 
376  // Open a stream and start writing scalar variable info.
377  std::ofstream scl_stream(scl_file.str().c_str());
378  scl_stream << "Per node scalar value\n";
379  scl_stream << "part\n";
380  scl_stream << std::setw(10) << 1 << "\n";
381  scl_stream << "coordinates\n";
382 
383  const MeshBase & the_mesh = MeshOutput<MeshBase>::mesh();
384  const unsigned int dim = the_mesh.mesh_dimension();
385  const System & system = _equation_systems.get_system(sys);
386  const DofMap & dof_map = system.get_dof_map();
387  int var = system.variable_number(var_name);
388 
389  std::vector<dof_id_type> dof_indices_scl;
390 
391  // Map from node id -> solution value. We end up just writing this
392  // map out in order, not sure what would happen if there were holes
393  // in the numbering...
394  typedef std::map<int, Real> map_local_soln;
395  typedef map_local_soln::iterator local_soln_iterator;
396  map_local_soln local_soln;
397 
398  std::vector<Number> elem_soln;
399  std::vector<Number> nodal_soln;
400 
401  // Loop over active local elements, construct the nodal solution, and write it to file.
402  for (const auto & elem : the_mesh.active_local_element_ptr_range())
403  {
404  const FEType & fe_type = system.variable_type(var);
405 
406  dof_map.dof_indices (elem, dof_indices_scl, var);
407 
408  elem_soln.resize(dof_indices_scl.size());
409 
410  for (std::size_t i = 0; i < dof_indices_scl.size(); i++)
411  elem_soln[i] = system.current_solution(dof_indices_scl[i]);
412 
413  FEInterface::nodal_soln (dim, fe_type, elem, elem_soln, nodal_soln);
414 
415  libmesh_assert_equal_to (nodal_soln.size(), elem->n_nodes());
416 
417 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
418  libmesh_error_msg("Complex-valued Ensight output not yet supported");
419 #endif
420 
421  for (unsigned int n=0; n<elem->n_nodes(); n++)
422  local_soln[elem->node_id(n)] = libmesh_real(nodal_soln[n]);
423  }
424 
425  {
426  local_soln_iterator it = local_soln.begin();
427  const local_soln_iterator it_end = local_soln.end();
428  for ( ; it != it_end; ++it)
429  scl_stream << std::setw(12)
430  << std::setprecision(5)
431  << std::scientific
432  << it->second
433  << "\n";
434  }
435 }
T libmesh_real(T a)
std::string _ensight_file_name
Definition: ensight_io.h:144
unsigned int dim
ImplicitSystem & sys
const MT & mesh() const
Definition: mesh_output.h:216
const EquationSystems & _equation_systems
Definition: ensight_io.h:152
const T_sys & get_system(const std::string &name) const
static void nodal_soln(const unsigned int dim, const FEType &fe_t, const Elem *elem, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln)
Build the nodal soln from the element soln.
Definition: fe_interface.C:526
std::vector< Real > _time_steps
Definition: ensight_io.h:145
void libMesh::EnsightIO::write_solution_ascii ( )
private

Definition at line 342 of file ensight_io.C.

References _system_vars_map, write_scalar_ascii(), and write_vector_ascii().

Referenced by write_ascii().

343 {
344  system_vars_map_t::iterator sys_it = _system_vars_map.begin();
345  const system_vars_map_t::iterator sys_end = _system_vars_map.end();
346 
347  for (; sys_it != sys_end; ++sys_it)
348  {
349  for (std::size_t i = 0; i < sys_it->second.EnsightScalars.size(); i++)
350  this->write_scalar_ascii(sys_it->first,
351  sys_it->second.EnsightScalars[i].scalar_name);
352 
353  for (std::size_t i = 0; i < sys_it->second.EnsightVectors.size(); i++)
354  this->write_vector_ascii(sys_it->first,
355  sys_it->second.EnsightVectors[i].components,
356  sys_it->second.EnsightVectors[i].description);
357  }
358 }
void write_scalar_ascii(const std::string &sys, const std::string &var)
Definition: ensight_io.C:361
system_vars_map_t _system_vars_map
Definition: ensight_io.h:149
void write_vector_ascii(const std::string &sys, const std::vector< std::string > &vec, const std::string &var_name)
Definition: ensight_io.C:438
void libMesh::EnsightIO::write_vector_ascii ( const std::string &  sys,
const std::vector< std::string > &  vec,
const std::string &  var_name 
)
private

Definition at line 438 of file ensight_io.C.

References _ensight_file_name, _equation_systems, _time_steps, libMesh::MeshBase::active_local_element_ptr_range(), libMesh::System::current_solution(), dim, libMesh::DofMap::dof_indices(), libMesh::System::get_dof_map(), libMesh::EquationSystems::get_system(), libMesh::libmesh_real(), libMesh::MeshOutput< MT >::mesh(), libMesh::MeshBase::mesh_dimension(), libMesh::FEInterface::nodal_soln(), libMesh::System::variable_number(), and libMesh::System::variable_type().

Referenced by write_solution_ascii().

441 {
442  // Construct vector variable filename
443  std::ostringstream vec_file;
444  vec_file << _ensight_file_name
445  << "_"
446  << var_name
447  << ".vec"
448  << std::setw(3)
449  << std::setprecision(0)
450  << std::setfill('0')
451  << std::right
452  << _time_steps.size()-1;
453 
454  // Open a stream and start writing vector variable info.
455  std::ofstream vec_stream(vec_file.str().c_str());
456  vec_stream << "Per vector per value\n";
457  vec_stream << "part\n";
458  vec_stream << std::setw(10) << 1 << "\n";
459  vec_stream << "coordinates\n";
460 
461  // Get a constant reference to the mesh object.
462  const MeshBase & the_mesh = MeshOutput<MeshBase>::mesh();
463 
464  // The dimension that we are running
465  const unsigned int dim = the_mesh.mesh_dimension();
466 
467  const System & system = _equation_systems.get_system(sys);
468 
469  const DofMap & dof_map = system.get_dof_map();
470 
471  const unsigned int u_var = system.variable_number(vec[0]);
472  const unsigned int v_var = system.variable_number(vec[1]);
473  const unsigned int w_var = (dim==3) ? system.variable_number(vec[2]) : 0;
474 
475  std::vector<dof_id_type> dof_indices_u;
476  std::vector<dof_id_type> dof_indices_v;
477  std::vector<dof_id_type> dof_indices_w;
478 
479  // Map from node id -> solution value. We end up just writing this
480  // map out in order, not sure what would happen if there were holes
481  // in the numbering...
482  typedef std::map<int,std::vector<Real>> map_local_soln;
483  typedef map_local_soln::iterator local_soln_iterator;
484  map_local_soln local_soln;
485 
486  // Now we will loop over all the elements in the mesh.
487  for (const auto & elem : the_mesh.active_local_element_ptr_range())
488  {
489  const FEType & fe_type = system.variable_type(u_var);
490 
491  dof_map.dof_indices (elem, dof_indices_u, u_var);
492  dof_map.dof_indices (elem, dof_indices_v, v_var);
493  if (dim==3)
494  dof_map.dof_indices (elem, dof_indices_w, w_var);
495 
496  std::vector<Number> elem_soln_u;
497  std::vector<Number> elem_soln_v;
498  std::vector<Number> elem_soln_w;
499 
500  std::vector<Number> nodal_soln_u;
501  std::vector<Number> nodal_soln_v;
502  std::vector<Number> nodal_soln_w;
503 
504  elem_soln_u.resize(dof_indices_u.size());
505  elem_soln_v.resize(dof_indices_v.size());
506  if (dim == 3)
507  elem_soln_w.resize(dof_indices_w.size());
508 
509  for (std::size_t i = 0; i < dof_indices_u.size(); i++)
510  {
511  elem_soln_u[i] = system.current_solution(dof_indices_u[i]);
512  elem_soln_v[i] = system.current_solution(dof_indices_v[i]);
513  if (dim==3)
514  elem_soln_w[i] = system.current_solution(dof_indices_w[i]);
515  }
516 
517  FEInterface::nodal_soln (dim, fe_type, elem, elem_soln_u, nodal_soln_u);
518  FEInterface::nodal_soln (dim, fe_type, elem, elem_soln_v, nodal_soln_v);
519  if (dim == 3)
520  FEInterface::nodal_soln (dim, fe_type, elem, elem_soln_w, nodal_soln_w);
521 
522  libmesh_assert_equal_to (nodal_soln_u.size(), elem->n_nodes());
523  libmesh_assert_equal_to (nodal_soln_v.size(), elem->n_nodes());
524 
525 #ifdef LIBMESH_ENABLE_COMPLEX
526  libmesh_error_msg("Complex-valued Ensight output not yet supported");
527 #endif
528 
529  for (unsigned int n=0; n<elem->n_nodes(); n++)
530  {
531  std::vector<Real> node_vec(3);
532  node_vec[0] = libmesh_real(nodal_soln_u[n]);
533  node_vec[1] = libmesh_real(nodal_soln_v[n]);
534  node_vec[2] = 0.0;
535  if (dim==3)
536  node_vec[2] = libmesh_real(nodal_soln_w[n]);
537  local_soln[elem->node_id(n)] = node_vec;
538  }
539  }
540 
541  {
542  local_soln_iterator it = local_soln.begin();
543  const local_soln_iterator it_end = local_soln.end();
544 
545  for (unsigned dir=0; dir<3; ++dir)
546  {
547  for (; it != it_end; ++it)
548  vec_stream << std::setw(12)
549  << std::scientific
550  << std::setprecision(5)
551  << it->second[dir]
552  << "\n";
553 
554  // Reset the iterator to the beginning of the map
555  it = local_soln.begin();
556  }
557  }
558 }
T libmesh_real(T a)
std::string _ensight_file_name
Definition: ensight_io.h:144
unsigned int dim
ImplicitSystem & sys
const MT & mesh() const
Definition: mesh_output.h:216
const EquationSystems & _equation_systems
Definition: ensight_io.h:152
const T_sys & get_system(const std::string &name) const
static void nodal_soln(const unsigned int dim, const FEType &fe_t, const Elem *elem, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln)
Build the nodal soln from the element soln.
Definition: fe_interface.C:526
std::vector< Real > _time_steps
Definition: ensight_io.h:145

Member Data Documentation

std::map< ElemType, std::string > libMesh::EnsightIO::_element_map = EnsightIO::build_element_map()
staticprivate

Definition at line 155 of file ensight_io.h.

Referenced by write_geometry_ascii().

std::string libMesh::EnsightIO::_ensight_file_name
private
const EquationSystems& libMesh::EnsightIO::_equation_systems
private
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 write().

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.

system_vars_map_t libMesh::EnsightIO::_system_vars_map
private

Definition at line 149 of file ensight_io.h.

Referenced by add_scalar(), add_vector(), write_case(), and write_solution_ascii().

std::vector<Real> libMesh::EnsightIO::_time_steps
private

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