#include <ensight_io.h>
Inherits MeshOutput< MeshBase >.
struct Scalars
struct SystemVars
struct Vectors
EnsightIO (const std::string &filename, const EquationSystems &eq)
~EnsightIO ()
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 add_scalar (const std::string &system, const std::string &scalar_description, const std::string &s)
virtual void write (const std::string &name)
void write (const double time=0)
bool & has_mesh_refinement ()
virtual void write_equation_systems (const std::string &, const EquationSystems &)
virtual void write_nodal_data (const std::string &, const std::vector< Number > &, const std::vector< std::string > &)
const MeshBase & mesh () const
typedef std::map< std::string, SystemVars > SystemsVarsMap
typedef std::map< std::string, SystemVars >::iterator SystemsVarsMapIterator
typedef std::pair< std::string, SystemVars > SystemsVarsValue
void write_ascii (const double 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 ()
void elem_type_to_string (ElemType, char *)
std::string _ensight_file_name
std::vector< double > _time_steps
SystemsVarsMap _systems_vars_map
const EquationSystems & _equation_systems
This class implements writing meshes and solutions in Ensight's Gold format. author Camata
Definition at line 42 of file ensight_io.h.
Definition at line 110 of file ensight_io.h.
Definition at line 111 of file ensight_io.h.
Definition at line 112 of file ensight_io.h.
Definition at line 37 of file ensight_io.C.
References _ensight_file_name, libMesh::n_processors(), and libMesh::processor_id().
:
MeshOutput<MeshBase> (eq.get_mesh()),
_equation_systems(eq)
{
if (libMesh::n_processors() == 1)
_ensight_file_name = filename;
else
{
std::stringstream tmp_file;
tmp_file << filename << '_rank' << libMesh::processor_id();
_ensight_file_name = tmp_file.str();
}
}
Definition at line 54 of file ensight_io.C.
{}
Definition at line 94 of file ensight_io.C.
References _equation_systems, _systems_vars_map, EnsightIO::Scalars::description, EquationSystems::get_system(), and EquationSystems::has_system().
{
libmesh_assert(_equation_systems.has_system(system_name));
libmesh_assert(_equation_systems.get_system(system_name).has_variable(s));
Scalars scl;
scl.description = scl_description;
scl.scalar_name = s;
_systems_vars_map[system_name].EnsightScalars.push_back(scl);
}
Definition at line 76 of file ensight_io.C.
References _equation_systems, _systems_vars_map, EnsightIO::Vectors::description, EquationSystems::get_system(), and EquationSystems::has_system().
{
libmesh_assert(_equation_systems.has_system(system_name));
libmesh_assert(_equation_systems.get_system(system_name).has_variable(u));
libmesh_assert(_equation_systems.get_system(system_name).has_variable(v));
libmesh_assert(_equation_systems.get_system(system_name).has_variable(w));
Vectors vec;
vec.description = vec_name;
vec.components.push_back(u);
vec.components.push_back(v);
vec.components.push_back(w);
_systems_vars_map[system_name].EnsightVectors.push_back(vec);
}
Definition at line 59 of file ensight_io.C.
References _equation_systems, _systems_vars_map, EnsightIO::Vectors::description, EquationSystems::get_system(), and EquationSystems::has_system().
{
libmesh_assert (_equation_systems.has_system(system_name));
libmesh_assert (_equation_systems.get_system(system_name).has_variable(u));
libmesh_assert (_equation_systems.get_system(system_name).has_variable(v));
Vectors vec;
vec.description = vec_description;
vec.components.push_back(u);
vec.components.push_back(v);
_systems_vars_map[system_name].EnsightVectors.push_back(vec);
}
Definition at line 537 of file ensight_io.C.
References libMeshEnums::EDGE2, libMeshEnums::EDGE3, libMeshEnums::HEX20, libMeshEnums::HEX27, libMeshEnums::HEX8, libMeshEnums::PYRAMID5, libMeshEnums::QUAD4, libMeshEnums::QUAD8, libMeshEnums::QUAD9, libMeshEnums::TET10, libMeshEnums::TET4, libMeshEnums::TRI3, and libMeshEnums::TRI6.
Referenced by write_geometry_ascii().
{
switch(type){
case EDGE2:
std::strcpy(buffer,'bar2');
break;
case EDGE3:
std::strcpy(buffer,'bar3');
break;
case QUAD4:
std::strcpy(buffer,'quad4');
break;
case QUAD8:
std::strcpy(buffer,'quad8');
break;
case QUAD9:
std::cout<<'QUAD9: element not supported!'<<std::endl;
libmesh_error();
break;
case TRI3:
std::strcpy(buffer,'tria3');
break;
case TRI6:
std::strcpy(buffer,'tria6');
break;
case TET4:
std::strcpy(buffer,'tetra4');
break;
case TET10:
std::strcpy(buffer,'tetra10');
break;
case HEX8:
std::strcpy(buffer,'hexa8');
break;
case HEX20:
std::strcpy(buffer,'hexa20');
break;
case HEX27:
std::cout<<'HEX27: element not supported!'<<std::endl;
libmesh_error();
break;
case PYRAMID5:
std::strcpy(buffer,'pyramid5');
break;
default:
break;
}
}
Referenced by PostscriptIO::write(), FroIO::write(), MEDITIO::write_ascii(), write_geometry_ascii(), write_scalar_ascii(), GnuPlotIO::write_solution(), DivaIO::write_stream(), and write_vector_ascii().
Implements MeshOutput< MeshBase >.
Definition at line 111 of file ensight_io.C.
References _ensight_file_name, and Quality::name().
{
_ensight_file_name = name;
this->write();
}
Definition at line 119 of file ensight_io.C.
References write_ascii(), and write_case().
{
this->write_ascii(time);
this->write_case();
}
Definition at line 128 of file ensight_io.C.
References _time_steps, write_geometry_ascii(), and write_solution_ascii().
Referenced by write().
{
_time_steps.push_back(time);
this->write_geometry_ascii();
this->write_solution_ascii();
}
Definition at line 260 of file ensight_io.C.
References _ensight_file_name, _systems_vars_map, _time_steps, EnsightIO::Vectors::description, EnsightIO::Scalars::description, and EnsightIO::Scalars::scalar_name.
Referenced by write().
{
std::stringstream case_file, geo_file, scl_file;
case_file << _ensight_file_name << '.case';
FILE* fout = fopen(case_file.str().c_str(),'w');
fprintf(fout,'FORMAT);
fprintf(fout,'type: ensight goldn');
fprintf(fout,'GEOMETRY);
geo_file << _ensight_file_name << '.geo';
fprintf(fout,'model: 1 %s***,geo_file.str().c_str());
SystemsVarsMapIterator sys = _systems_vars_map.begin();
const SystemsVarsMapIterator sys_end = _systems_vars_map.end();
// Write Variable per node section
if ( sys != sys_end )
fprintf(fout,'nVARIABLE);
for (; sys != sys_end; ++sys)
{
SystemsVarsValue value = *sys;
for (unsigned int i=0; i < value.second.EnsightScalars.size(); i++)
{
std::stringstream scl_file;
Scalars scalar = value.second.EnsightScalars[i];
scl_file << _ensight_file_name
<< '_' << scalar.scalar_name
<< '.scl';
fprintf(fout,'scalar per node: 1 %s %s***,scalar.description.c_str(), scl_file.str().c_str());
}
for (unsigned int i=0; i < value.second.EnsightVectors.size(); i++)
{
std::stringstream vec_file;
Vectors vec = value.second.EnsightVectors[i];
vec_file<<_ensight_file_name<<'_'<<vec.description<<'.vec';
fprintf(fout,'vector per node: 1 %s %s***,vec.description.c_str(), vec_file.str().c_str());
}
// Write time step section
if( _time_steps.size() != 0)
{
fprintf(fout,'nTIME);
fprintf(fout,'time set: 1);
fprintf(fout,'number of steps: %10d, static_cast<int>(_time_steps.size()));
fprintf(fout,'filename start number: %10d, 0);
fprintf(fout,'filename increment: %10d, 1);
fprintf(fout,'time values:);
for (unsigned int i = 0; i < _time_steps.size(); i++)
fprintf(fout,'%12.5e, _time_steps[i]);
}
}
fclose(fout);
}
Definition at line 138 of file ensight_io.C.
References _ensight_file_name, _time_steps, MeshBase::active_local_elements_begin(), MeshBase::active_local_elements_end(), elem_type_to_string(), libMeshEnums::HEX27, MeshOutput< MeshBase >::mesh(), Elem::n_nodes(), Elem::node(), Elem::point(), libMeshEnums::QUAD9, and Elem::type().
Referenced by write_ascii().
{
OStringStream file;
file << _ensight_file_name << '.geo';
OSSRealzeroright(file,3,0,_time_steps.size()-1);
FILE* fout = fopen(file.str().c_str(),'w');
char buffer[80];
fprintf(fout,'EnSight Gold Geometry File Format);
fprintf(fout,'Generated by );
fprintf(fout,'node id off);
fprintf(fout,'element id given);
fprintf(fout,'part);
fprintf(fout,'%10d,1);
fprintf(fout,'uns-elements);
fprintf(fout,'coordinates);
// mapping between nodal index and your coordinates
std::map<int, Point> mesh_nodes_map;
typedef std::map <int, Point>::iterator mesh_nodes_iterator;
typedef std::pair<int, Point> mesh_node_value;
// Mapping between global and local indices
std::map <int, int> ensight_node_index;
// Grouping elements of the same type
std::map<ElemType, std::vector<const Elem*> > ensight_parts_map;
typedef std::map<ElemType, std::vector<const Elem*> >::iterator ensight_parts_iterator;
typedef std::pair<ElemType, std::vector<const Elem*> > ensight_parts_value;
const MeshBase& mesh = MeshOutput<MeshBase>::mesh();
MeshBase::const_element_iterator el = mesh.active_local_elements_begin();
const MeshBase::const_element_iterator end_el = mesh.active_local_elements_end();
for ( ; el != end_el ; ++el)
{
const Elem* elem = *el;
ensight_parts_map[elem->type()].push_back(elem);
for (unsigned int i = 0; i < elem->n_nodes(); i++)
mesh_nodes_map[elem->node(i)] = elem->point(i);
}
// Write number of local points
fprintf(fout,'%10d,static_cast<int>(mesh_nodes_map.size()));
mesh_nodes_iterator no_it = mesh_nodes_map.begin();
const mesh_nodes_iterator no_end_it = mesh_nodes_map.end();
// write x
for(int i = 1; no_it != no_end_it; ++no_it, i++)
{
const mesh_node_value pn = *no_it;
fprintf(fout,'%12.5e,pn.second(0));
ensight_node_index[pn.first] = i;
}
// write y
no_it = mesh_nodes_map.begin();
for(; no_it != no_end_it; ++no_it)
{
const mesh_node_value pn = *no_it;
fprintf(fout,'%12.5e,pn.second(1));
}
// write z
no_it = mesh_nodes_map.begin();
for(; no_it != no_end_it; ++no_it)
{
const mesh_node_value pn = *no_it;
fprintf(fout,'%12.5e,pn.second(2));
}
ensight_parts_iterator parts_it = ensight_parts_map.begin();
const ensight_parts_iterator end_parts_it = ensight_parts_map.end();
// Write parts
for (; parts_it != end_parts_it; ++parts_it)
{
ensight_parts_value kvp = *parts_it;
// Write element type
elem_type_to_string(kvp.first,buffer);
fprintf(fout,'s, buffer);
std::vector<const Elem*> elem_ref = kvp.second;
// Write number of element
fprintf(fout,'%10d,static_cast<int>(elem_ref.size()));
// Write element id
for (unsigned int i = 0; i < elem_ref.size(); i++)
fprintf(fout,'%10d,elem_ref[i]->id());
// Write connectivity
for (unsigned int i = 0; i < elem_ref.size(); i++)
{
for (unsigned int j = 0; j < elem_ref[i]->n_nodes(); j++) {
// tests!
if(kvp.first == QUAD9 && i==4)
continue;
// tests!
if(kvp.first == HEX27 && (i==4 || i ==10 || i == 12 ||
i == 13 || i ==14 || i == 16 || i == 22))
continue;
fprintf(fout,'%10d',ensight_node_index[elem_ref[i]->node(j)]);
}
fprintf(fout,');
}
}
fclose(fout);
}
Reimplemented in ExodusII_IO, GmshIO, GMVIO, GnuPlotIO, MEDITIO, and TecplotIO.
Definition at line 91 of file mesh_output.h.
{ libmesh_error(); }
Definition at line 346 of file ensight_io.C.
References _ensight_file_name, _equation_systems, _time_steps, MeshBase::active_local_elements_begin(), MeshBase::active_local_elements_end(), System::current_solution(), DofMap::dof_indices(), System::get_dof_map(), EquationSystems::get_system(), libmesh_real(), MeshOutput< MeshBase >::mesh(), MeshBase::mesh_dimension(), Elem::n_nodes(), FEInterface::nodal_soln(), Elem::node(), System::variable_number(), and System::variable_type().
Referenced by write_solution_ascii().
{
OStringStream scl_file;
scl_file << _ensight_file_name << '_' << var_name << '.scl';
OSSRealzeroright(scl_file,3,0,_time_steps.size()-1);
FILE * fout = fopen(scl_file.str().c_str(),'w');
fprintf(fout,'Per node scalar value);
fprintf(fout,'part);
fprintf(fout,'%10d,1);
fprintf(fout,'coordinates);
const MeshBase& mesh = MeshOutput<MeshBase>::mesh();
const unsigned int dim = mesh.mesh_dimension();
const System &system = _equation_systems.get_system(sys);
const DofMap& dof_map = system.get_dof_map();
int var = system.variable_number(var_name);
std::vector<unsigned int> dof_indices;
std::vector<unsigned int> dof_indices_scl;
// Now we will loop over all the elements in the mesh.
MeshBase::const_element_iterator el = mesh.active_local_elements_begin();
const MeshBase::const_element_iterator end_el = mesh.active_local_elements_end();
typedef std::map<int,Real> map_local_soln;
typedef map_local_soln::iterator local_soln_iterator;
map_local_soln local_soln;
std::vector<Number> elem_soln;
std::vector<Number> nodal_soln;
for ( ; el != end_el ; ++el){
const Elem* elem = *el;
const FEType& fe_type = system.variable_type(var);
dof_map.dof_indices (elem, dof_indices);
dof_map.dof_indices (elem, dof_indices_scl, var);
elem_soln.resize(dof_indices_scl.size());
for (unsigned int i = 0; i < dof_indices_scl.size(); i++)
elem_soln[i] = system.current_solution(dof_indices_scl[i]);
FEInterface::nodal_soln (dim,fe_type, elem, elem_soln, nodal_soln);
libmesh_assert (nodal_soln.size() == elem->n_nodes());
#ifdef LIBMESH_ENABLE_COMPLEX
std::cerr << 'Complex-valued Ensight output not yet supported' << std::endl;
libmesh_not_implemented()
#endif
for (unsigned int n=0; n<elem->n_nodes(); n++)
local_soln[elem->node(n)] = libmesh_real(nodal_soln[n]);
}
local_soln_iterator sol = local_soln.begin();
const local_soln_iterator sol_end = local_soln.end();
for(; sol != sol_end; ++sol)
fprintf(fout,'%12.5e,(*sol).second);
fclose(fout);
}
Definition at line 324 of file ensight_io.C.
References _systems_vars_map, write_scalar_ascii(), and write_vector_ascii().
Referenced by write_ascii().
{
SystemsVarsMapIterator sys = _systems_vars_map.begin();
const SystemsVarsMapIterator sys_end = _systems_vars_map.end();
for (; sys != sys_end; ++sys)
{
SystemsVarsValue value = *sys;
for (unsigned int i = 0; i < value.second.EnsightScalars.size(); i++)
this->write_scalar_ascii(value.first,
value.second.EnsightScalars[i].scalar_name);
for (unsigned int i = 0; i < value.second.EnsightVectors.size(); i++)
this->write_vector_ascii(value.first,
value.second.EnsightVectors[i].components,
value.second.EnsightVectors[i].description);
}
}
Definition at line 425 of file ensight_io.C.
References _ensight_file_name, _equation_systems, _time_steps, MeshBase::active_local_elements_begin(), MeshBase::active_local_elements_end(), System::current_solution(), DofMap::dof_indices(), System::get_dof_map(), EquationSystems::get_system(), libmesh_real(), MeshOutput< MeshBase >::mesh(), MeshBase::mesh_dimension(), Elem::n_nodes(), FEInterface::nodal_soln(), Elem::node(), System::variable_number(), and System::variable_type().
Referenced by write_solution_ascii().
{
OStringStream vec_file;
vec_file<<_ensight_file_name<<'_'<<var_name<<'.vec';
OSSRealzeroright(vec_file,3,0,_time_steps.size()-1);
FILE * fout = fopen(vec_file.str().c_str(),'w');
fprintf(fout,'Per vector per value);
fprintf(fout,'part);
fprintf(fout,'%10d,1);
fprintf(fout,'coordinates);
// Get a constant reference to the mesh object.
const MeshBase& mesh = MeshOutput<MeshBase>::mesh();
// The dimension that we are running
const unsigned int dim = mesh.mesh_dimension();
const System &system = _equation_systems.get_system(sys);
const DofMap& dof_map = system.get_dof_map();
const unsigned int u_var = system.variable_number(vec[0]);
const unsigned int v_var = system.variable_number(vec[1]);
const unsigned int w_var = (dim==3) ? system.variable_number(vec[2]) : 0;
std::vector<unsigned int> dof_indices;
std::vector<unsigned int> dof_indices_u;
std::vector<unsigned int> dof_indices_v;
std::vector<unsigned int> dof_indices_w;
// Now we will loop over all the elements in the mesh.
MeshBase::const_element_iterator el = mesh.active_local_elements_begin();
const MeshBase::const_element_iterator end_el = mesh.active_local_elements_end();
typedef std::map<int,std::vector<Real> > map_local_soln;
typedef map_local_soln::iterator local_soln_iterator;
map_local_soln local_soln;
for ( ; el != end_el ; ++el){
const Elem* elem = *el;
const FEType& fe_type = system.variable_type(u_var);
dof_map.dof_indices (elem, dof_indices);
dof_map.dof_indices (elem, dof_indices_u,u_var);
dof_map.dof_indices (elem, dof_indices_v,v_var);
if(dim==3) dof_map.dof_indices (elem, dof_indices,w_var);
std::vector<Number> elem_soln_u;
std::vector<Number> elem_soln_v;
std::vector<Number> elem_soln_w;
std::vector<Number> nodal_soln_u;
std::vector<Number> nodal_soln_v;
std::vector<Number> nodal_soln_w;
elem_soln_u.resize(dof_indices_u.size());
elem_soln_v.resize(dof_indices_v.size());
if(dim == 3) elem_soln_w.resize(dof_indices_w.size());
for (unsigned int i = 0; i < dof_indices_u.size(); i++)
{
elem_soln_u[i] = system.current_solution(dof_indices_u[i]);
elem_soln_v[i] = system.current_solution(dof_indices_v[i]);
if(dim==3) elem_soln_w[i] = system.current_solution(dof_indices_w[i]);
}
FEInterface::nodal_soln (dim,fe_type,elem,elem_soln_u,nodal_soln_u);
FEInterface::nodal_soln (dim,fe_type,elem,elem_soln_v,nodal_soln_v);
if(dim == 3) FEInterface::nodal_soln (dim,fe_type,elem,elem_soln_w,nodal_soln_w);
libmesh_assert (nodal_soln_u.size() == elem->n_nodes());
libmesh_assert (nodal_soln_v.size() == elem->n_nodes());
#ifdef LIBMESH_ENABLE_COMPLEX
std::cerr << 'Complex-valued Ensight output not yet supported' << std::endl;
libmesh_not_implemented()
#endif
for (unsigned int n=0; n<elem->n_nodes(); n++)
{
std::vector<Real> vec(3);
vec[0]= libmesh_real(nodal_soln_u[n]);
vec[1]= libmesh_real(nodal_soln_v[n]);
vec[2]=0.0;
if(dim==3) vec[2]= libmesh_real(nodal_soln_w[n]);
local_soln[elem->node(n)] = vec;
}
}
local_soln_iterator sol = local_soln.begin();
const local_soln_iterator sol_end = local_soln.end();
for(; sol != sol_end; ++sol)
fprintf(fout,'%12.5e,(*sol).second[0]);
sol = local_soln.begin();
for(; sol != sol_end; ++sol)
fprintf(fout,'%12.5e,(*sol).second[1]);
sol = local_soln.begin();
for(; sol != sol_end; ++sol)
fprintf(fout,'%12.5e,(*sol).second[2]);
fclose(fout);
}
Definition at line 127 of file ensight_io.h.
Referenced by EnsightIO(), write(), write_case(), write_geometry_ascii(), write_scalar_ascii(), and write_vector_ascii().
Definition at line 130 of file ensight_io.h.
Referenced by add_scalar(), add_vector(), write_scalar_ascii(), and write_vector_ascii().
Definition at line 129 of file ensight_io.h.
Referenced by add_scalar(), add_vector(), write_case(), and write_solution_ascii().
Definition at line 128 of file ensight_io.h.
Referenced by write_ascii(), write_case(), write_geometry_ascii(), write_scalar_ascii(), and write_vector_ascii().
Generated automatically by Doxygen for libMesh from the source code.