libMesh::FroIO Class Reference

#include <fro_io.h>

Inheritance diagram for libMesh::FroIO:

Public Member Functions

 FroIO (const MeshBase &)
 
virtual void write (const std::string &)
 
virtual void write_equation_systems (const std::string &, const EquationSystems &, const std::set< std::string > *system_names=NULL)
 
virtual void write_nodal_data (const std::string &, const std::vector< Number > &, const std::vector< std::string > &)
 
unsigned int & ascii_precision ()
 

Protected Member Functions

const MeshBasemesh () const
 

Protected Attributes

const bool _is_parallel_format
 

Detailed Description

This class implements writing meshes in the .fro format used by the MIT ACDL. Valid only for triangular meshes.

Author
Benjamin S. Kirk, 2007

Definition at line 46 of file fro_io.h.

Constructor & Destructor Documentation

libMesh::FroIO::FroIO ( const MeshBase mesh_in)
inlineexplicit

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

Definition at line 72 of file fro_io.h.

72  :
73  MeshOutput<MeshBase> (mesh_in)
74 {
75 }

Member Function Documentation

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

void libMesh::FroIO::write ( const std::string &  fname)
virtual

This method implements writing a mesh to a specified file.

Implements libMesh::MeshOutput< MeshBase >.

Definition at line 41 of file fro_io.C.

References libMesh::MeshOutput< MeshBase >::_is_parallel_format, libMesh::MeshBase::active_elements_begin(), libMesh::MeshBase::active_elements_end(), bc_id, libMesh::MeshBase::boundary_info, libMesh::Elem::build_side(), libMesh::MeshBase::elem(), end, libMesh::err, libMesh::libmesh_assert(), libMesh::MeshOutput< MT >::mesh(), libMesh::MeshOutput< MeshBase >::mesh(), libMesh::MeshBase::n_elem(), libMesh::MeshBase::n_nodes(), libMesh::MeshBase::point(), libMesh::processor_id(), side, and libMeshEnums::TRI3.

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

42 {
43  // We may need to gather a ParallelMesh to output it, making that
44  // const qualifier in our constructor a dirty lie
45  MeshSerializer serialize(const_cast<MeshBase&>(this->mesh()), !_is_parallel_format);
46 
47  if (this->mesh().processor_id() == 0)
48  {
49  // Open the output file stream
50  std::ofstream out_stream (fname.c_str());
51  libmesh_assert (out_stream.good());
52 
53  // Make sure it opened correctly
54  if (!out_stream.good())
55  libmesh_file_error(fname.c_str());
56 
57  // Get a reference to the mesh
58  const MeshBase& the_mesh = MeshOutput<MeshBase>::mesh();
59 
60  // Write the header
61  out_stream << the_mesh.n_elem() << " "
62  << the_mesh.n_nodes() << " "
63  << "0 0 "
64  << the_mesh.boundary_info->n_boundary_ids() << " 1\n";
65 
66  // Write the nodes -- 1-based!
67  for (unsigned int n=0; n<the_mesh.n_nodes(); n++)
68  out_stream << n+1 << " \t"
69  << std::scientific
70  << std::setprecision(12)
71  << the_mesh.point(n)(0) << " \t"
72  << the_mesh.point(n)(1) << " \t"
73  << 0. << '\n';
74 
75  // Write the elements -- 1-based!
76  MeshBase::const_element_iterator it = the_mesh.active_elements_begin();
77  const MeshBase::const_element_iterator end = the_mesh.active_elements_end();
78 
79  for (unsigned int e=0 ; it != end; ++it)
80  {
81  // .fro likes TRI3's
82  if ((*it)->type() != TRI3)
83  {
84  libMesh::err << "ERROR: .fro format only valid for triangles!\n"
85  << " writing of " << fname << " aborted.\n"
86  << std::endl;
87  libmesh_error();
88  }
89 
90  out_stream << ++e << " \t";
91 
92  for (unsigned int n=0; n<(*it)->n_nodes(); n++)
93  out_stream << (*it)->node(n)+1 << " \t";
94 
95 // // LHS -> RHS Mapping, for inverted triangles
96 // out_stream << (*it)->node(0)+1 << " \t";
97 // out_stream << (*it)->node(2)+1 << " \t";
98 // out_stream << (*it)->node(1)+1 << " \t";
99 
100  out_stream << "1\n";
101  }
102 
103  // Write BCs.
104  {
105  const std::set<boundary_id_type>& bc_ids =
106  the_mesh.boundary_info->get_boundary_ids();
107 
108  std::vector<dof_id_type> el;
109  std::vector<unsigned short int> sl;
110  std::vector<boundary_id_type> il;
111 
112  the_mesh.boundary_info->build_side_list (el, sl, il);
113 
114 
115  // Map the boundary ids into [1,n_bc_ids],
116  // treat them one at a time.
118  for (std::set<boundary_id_type>::const_iterator id = bc_ids.begin();
119  id != bc_ids.end(); ++id)
120  {
121  std::deque<dof_id_type> node_list;
122 
123  std::map<dof_id_type, dof_id_type>
124  forward_edges, backward_edges;
125 
126  // Get all sides on this element with the relevant BC id.
127  for (std::size_t e=0; e<el.size(); e++)
128  if (il[e] == *id)
129  {
130  // need to build up node_list as a sorted array of edge nodes...
131  // for the following:
132  // a---b---c---d---e
133  // node_list [ a b c d e];
134  //
135  // the issue is just how to get this out of the elem/side based data structure.
136  // the approach is to build up 'chain links' like this:
137  // a---b b---c c---d d---e
138  // and piece them together.
139  //
140  // so, for an arbitray edge n0---n1, we build the
141  // "forward_edges" map n0-->n1
142  // "backward_edges" map n1-->n0
143  // and then start with one chain link, and add on...
144  //
145  AutoPtr<Elem> side = the_mesh.elem(el[e])->build_side(sl[e]);
146 
147  const dof_id_type
148  n0 = side->node(0),
149  n1 = side->node(1);
150 
151  // insert into forward-edge set
152  forward_edges.insert (std::make_pair(n0, n1));
153 
154  // insert into backward-edge set
155  backward_edges.insert (std::make_pair(n1, n0));
156 
157  // go ahead and add one edge to the list -- this will give us the beginning of a
158  // chain to work from!
159  if (node_list.empty())
160  {
161  node_list.push_front(n0);
162  node_list.push_back (n1);
163  }
164  }
165 
166  // we now have the node_list with one edge, the forward_edges, and the backward_edges
167  // the node_list will be filled when (node_list.size() == (n_edges+1))
168  // until that is the case simply add on to the beginning and end of the node_list,
169  // building up a chain of ordered nodes...
170  const std::size_t n_edges = forward_edges.size();
171 
172  while (node_list.size() != (n_edges+1))
173  {
174  const dof_id_type
175  front_node = node_list.front(),
176  back_node = node_list.back();
177 
178  // look for front_pair in the backward_edges list
179  {
180  std::map<dof_id_type, dof_id_type>::iterator
181  pos = backward_edges.find(front_node);
182 
183  if (pos != backward_edges.end())
184  {
185  node_list.push_front(pos->second);
186 
187  backward_edges.erase(pos);
188  }
189  }
190 
191  // look for back_pair in the forward_edges list
192  {
193  std::map<dof_id_type, dof_id_type>::iterator
194  pos = forward_edges.find(back_node);
195 
196  if (pos != forward_edges.end())
197  {
198  node_list.push_back(pos->second);
199 
200  forward_edges.erase(pos);
201  }
202  }
203 
204 // libMesh::out << "node_list.size()=" << node_list.size()
205 // << ", n_edges+1=" << n_edges+1 << std::endl;
206  }
207 
208 
209  out_stream << ++bc_id << " " << node_list.size() << '\n';
210 
211  std::deque<dof_id_type>::iterator pos = node_list.begin();
212  for ( ; pos != node_list.end(); ++pos)
213  out_stream << *pos+1 << " \t0\n";
214  }
215  }
216  }
217 }
virtual void libMesh::MeshOutput< MeshBase >::write_equation_systems ( const std::string &  ,
const EquationSystems ,
const std::set< std::string > *  system_names = NULL 
)
virtualinherited

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

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

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::GMVIO, libMesh::Nemesis_IO, libMesh::GmshIO, libMesh::VTKIO, libMesh::UCDIO, libMesh::MEDITIO, libMesh::GnuPlotIO, and libMesh::TecplotIO.

Definition at line 98 of file mesh_output.h.

101  { libmesh_error(); }

Member Data Documentation

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

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


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

Site Created By: libMesh Developers
Last modified: February 07 2014 16:58:01 UTC

Hosted By:
SourceForge.net Logo