libMesh::EnsightIO Class Reference

#include <ensight_io.h>

Inheritance diagram for libMesh::EnsightIO:

Classes

struct  Scalars
 
struct  SystemVars
 
struct  Vectors
 

Public Member Functions

 EnsightIO (const std::string &filename, const EquationSystems &eq)
 
 ~EnsightIO ()
 
void add_scalar (const std::string &system, const std::string &scalar_description, const std::string &s)
 
void add_vector (const std::string &system, const std::string &vec_description, const std::string &u, const std::string &v)
 
void add_vector (const std::string &system, const std::string &vec_description, const std::string &u, const std::string &v, const std::string &w)
 
void write (Real time=0)
 
virtual void write (const std::string &name) libmesh_override
 
virtual void write_equation_systems (const std::string &, const EquationSystems &, const std::set< std::string > *system_names=libmesh_nullptr)
 
virtual void write_nodal_data (const std::string &, const std::vector< Number > &, const std::vector< std::string > &)
 
virtual void write_nodal_data (const std::string &, const NumericVector< Number > &, const std::vector< std::string > &)
 
unsigned int & ascii_precision ()
 

Protected Member Functions

const MeshBasemesh () const
 

Protected Attributes

const bool _is_parallel_format
 
const bool _serial_only_needed_on_proc_0
 

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

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 }
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)
Definition: elem_quality.C:39
std::string _ensight_file_name
Definition: ensight_io.h:144
const MeshBase & mesh() const
void write(Real time=0)
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:348
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 293 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().

294 {
295  std::ostringstream case_file;
296  case_file << _ensight_file_name << ".case";
297 
298  // Open a stream for writing the case file.
299  std::ofstream case_stream(case_file.str().c_str());
300 
301  case_stream << "FORMAT\n";
302  case_stream << "type: ensight gold\n\n";
303  case_stream << "GEOMETRY\n";
304  case_stream << "model: 1 " << _ensight_file_name << ".geo" << "***\n";
305 
306  system_vars_map_t::iterator sys_it = _system_vars_map.begin();
307  const system_vars_map_t::iterator sys_end = _system_vars_map.end();
308 
309  // Write Variable per node section
310  if (sys_it != sys_end)
311  case_stream << "\n\nVARIABLE\n";
312 
313  for (; sys_it != sys_end; ++sys_it)
314  {
315  for (std::size_t i=0; i < sys_it->second.EnsightScalars.size(); i++)
316  {
317  Scalars scalar = sys_it->second.EnsightScalars[i];
318  case_stream << "scalar per node: 1 "
319  << scalar.description << " "
320  << _ensight_file_name << "_" << scalar.scalar_name << ".scl***\n";
321  }
322 
323  for (std::size_t i=0; i < sys_it->second.EnsightVectors.size(); i++)
324  {
325  Vectors vec = sys_it->second.EnsightVectors[i];
326  case_stream << "vector per node: 1 "
327  << vec.description << " "
328  << _ensight_file_name << "_" << vec.description << ".vec***\n";
329  }
330 
331  // Write time step section
332  if (_time_steps.size() != 0)
333  {
334  case_stream << "\n\nTIME\n";
335  case_stream << "time set: 1\n";
336  case_stream << "number of steps: " << std::setw(10) << _time_steps.size() << "\n";
337  case_stream << "filename start number: " << std::setw(10) << 0 << "\n";
338  case_stream << "filename increment: " << std::setw(10) << 1 << "\n";
339  case_stream << "time values:\n";
340  for (std::size_t i = 0; i < _time_steps.size(); i++)
341  case_stream << std::setw(12) << std::setprecision(5) << std::scientific << _time_steps[i] << "\n";
342  }
343  }
344 }
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_elements_begin(), libMesh::MeshBase::active_local_elements_end(), libMesh::HEX27, libMesh::MeshOutput< MT >::mesh(), libMesh::Elem::n_nodes(), libMesh::Elem::node_id(), libMesh::Elem::point(), libMesh::QUAD9, libMesh::TypeVector< T >::size(), and libMesh::Elem::type().

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  {
202  MeshBase::const_element_iterator el = the_mesh.active_local_elements_begin();
203  const MeshBase::const_element_iterator end_el = the_mesh.active_local_elements_end();
204 
205  for ( ; el != end_el ; ++el)
206  {
207  const Elem * elem = *el;
208  ensight_parts_map[elem->type()].push_back(elem);
209 
210  for (unsigned int i = 0; i < elem->n_nodes(); i++)
211  mesh_nodes_map[elem->node_id(i)] = elem->point(i);
212  }
213  }
214 
215  // Write number of local points
216  mesh_stream << std::setw(10) << mesh_nodes_map.size() << "\n";
217 
218  // write x, y, and z node positions, build mapping between
219  // ensight and libmesh node numbers.
220  std::map <int, int> ensight_node_index;
221  {
222  mesh_nodes_iterator no_it = mesh_nodes_map.begin();
223  const mesh_nodes_iterator no_end_it = mesh_nodes_map.end();
224 
225  for (unsigned direction=0; direction<3; ++direction)
226  {
227  for (int i = 1; no_it != no_end_it; ++no_it, i++)
228  {
229  mesh_stream << std::setw(12)
230  << std::setprecision(5)
231  << std::scientific
232  << no_it->second(direction)
233  << "\n";
234  ensight_node_index[no_it->first] = i;
235  }
236 
237  // Reset iterator to the beginning of the map
238  no_it = mesh_nodes_map.begin();
239  }
240  }
241 
242  // Write parts
243  {
244  ensight_parts_iterator parts_it = ensight_parts_map.begin();
245  const ensight_parts_iterator end_parts_it = ensight_parts_map.end();
246 
247  for (; parts_it != end_parts_it; ++parts_it)
248  {
249  // Look up this ElemType in the map, error if not present.
250  std::map<ElemType, std::string>::iterator name_it = _element_map.find(parts_it->first);
251  if (name_it == _element_map.end())
252  libmesh_error_msg("Error: Unsupported ElemType " << parts_it->first << " for EnsightIO.");
253 
254  // Write element type
255  mesh_stream << "\n" << name_it->second << "\n";
256 
257  std::vector<const Elem *> elem_ref = parts_it->second;
258 
259  // Write number of element
260  mesh_stream << std::setw(10) << elem_ref.size() << "\n";
261 
262  // Write element id
263  for (std::size_t i = 0; i < elem_ref.size(); i++)
264  mesh_stream << std::setw(10) << elem_ref[i]->id() << "\n";
265 
266  // Write connectivity
267  for (std::size_t i = 0; i < elem_ref.size(); i++)
268  {
269  for (unsigned int j = 0; j < elem_ref[i]->n_nodes(); j++)
270  {
271  // tests!
272  if (parts_it->first == QUAD9 && i==4)
273  continue;
274 
275  // tests!
276  if (parts_it->first == HEX27 &&
277  (i==4 || i ==10 || i == 12 ||
278  i == 13 || i ==14 || i == 16 || i == 22))
279  continue;
280 
281  mesh_stream << std::setw(10) << ensight_node_index[elem_ref[i]->node_id(j)];
282  }
283  mesh_stream << "\n";
284  }
285  }
286  }
287 }
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 > &   
)
inlinevirtualinherited

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 367 of file ensight_io.C.

References _ensight_file_name, _equation_systems, _time_steps, libMesh::MeshBase::active_local_elements_begin(), libMesh::MeshBase::active_local_elements_end(), libMesh::System::current_solution(), 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::Elem::n_nodes(), libMesh::FEInterface::nodal_soln(), libMesh::Elem::node_id(), libMesh::System::variable_number(), and libMesh::System::variable_type().

Referenced by write_solution_ascii().

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

Definition at line 348 of file ensight_io.C.

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

Referenced by write_ascii().

349 {
350  system_vars_map_t::iterator sys_it = _system_vars_map.begin();
351  const system_vars_map_t::iterator sys_end = _system_vars_map.end();
352 
353  for (; sys_it != sys_end; ++sys_it)
354  {
355  for (std::size_t i = 0; i < sys_it->second.EnsightScalars.size(); i++)
356  this->write_scalar_ascii(sys_it->first,
357  sys_it->second.EnsightScalars[i].scalar_name);
358 
359  for (std::size_t i = 0; i < sys_it->second.EnsightVectors.size(); i++)
360  this->write_vector_ascii(sys_it->first,
361  sys_it->second.EnsightVectors[i].components,
362  sys_it->second.EnsightVectors[i].description);
363  }
364 }
void write_scalar_ascii(const std::string &sys, const std::string &var)
Definition: ensight_io.C:367
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:449
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 449 of file ensight_io.C.

References _ensight_file_name, _equation_systems, _time_steps, libMesh::MeshBase::active_local_elements_begin(), libMesh::MeshBase::active_local_elements_end(), libMesh::System::current_solution(), 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::Elem::n_nodes(), libMesh::FEInterface::nodal_soln(), libMesh::Elem::node_id(), libMesh::System::variable_number(), and libMesh::System::variable_type().

Referenced by write_solution_ascii().

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