*/
namespace GridOutFlags
{
+ /**
+ * Flags for grid output in OpenDX format.
+ */
+ struct DX
+ {
+ /**
+ * Write faces instead of
+ * cells. Defaults to false.
+ */
+ bool write_faces;
+ /**
+ * Write all faces, including
+ * interior faces. Default is to
+ * write boundary faces only.
+ */
+ bool write_all_faces;
+
+ /**
+ * Constructor.
+ */
+ DX (bool write_faces = false,
+ bool write_all_faces = false);
+ };
+
/**
* Flags describing the details of
* output in UCD format.
* Declaration of a name for each of the
* different output formats.
*/
- enum OutputFormat { gnuplot, eps, ucd };
+ enum OutputFormat { dx, gnuplot, eps, ucd };
+
+ /**
+ * Write triangulation in OpenDX format.
+ *
+ * Cells or faces are written
+ * together with their level and
+ * their material id or boundary
+ * indicator, resp.
+ */
+ template <int dim>
+ void write_dx (const Triangulation<dim> &tria,
+ std::ostream &out);
/**
* Write the triangulation in the
const Mapping<dim> *mapping=0);
/**
- * Set the flags to be used for output
- * in UCD format.
+ * Set flags for DX output
*/
- void set_flags (const GridOutFlags::Ucd &ucd_flags);
+ void set_flags (const GridOutFlags::DX &flags);
+
+ /**
+ * Set flags for UCD output
+ */
+ void set_flags (const GridOutFlags::Ucd &flags);
/**
- * Set the flags to be used for output
- * in GNUPLOT format.
+ * Set flags for GNUPLOT output
*/
- void set_flags (const GridOutFlags::Gnuplot &gnuplot_flags);
+ void set_flags (const GridOutFlags::Gnuplot &flags);
/**
- * Set the flags to be used for output
- * in 1d EPS output.
+ * Set flags for EPS output of a
+ * one-dimensional triangulation
*/
- void set_flags (const GridOutFlags::Eps<1> &eps_flags);
+ void set_flags (const GridOutFlags::Eps<1> &flags);
- /**
- * Set the flags to be used for output
- * in 2d EPS output.
+ /**
+ * Set flags for EPS output of a
+ * two-dimensional triangulation
*/
- void set_flags (const GridOutFlags::Eps<2> &eps_flags);
+ void set_flags (const GridOutFlags::Eps<2> &flags);
- /**
- * Set the flags to be used for output
- * in 3d EPS output.
+ /**
+ * Set flags for EPS output of a
+ * three-dimensional triangulation
*/
- void set_flags (const GridOutFlags::Eps<3> &eps_flags);
+ void set_flags (const GridOutFlags::Eps<3> &flags);
/**
* Provide a function which tells us which
* usually has. At present the following
* formats are defined:
* @begin{itemize}
+ * @item @p{OpenDX}: @p{.dx}
* @item @p{gnuplot}: @p{.gnuplot}
* @item @p{ucd}: @p{.inp}
* @item @p{eps}: @p{.eps}.
DeclException0 (ExcInvalidState);
private:
-
+ /**
+ * Flags for OpenDX output.
+ */
+ GridOutFlags::DX dx_flags;
+
/**
- * Flags to be used upon output of UCD
- * data. Can be changed by using the
+ * Flags for UCD output. Can be
+ * changed by using the
* @p{set_flags} function.
*/
GridOutFlags::Ucd ucd_flags;
namespace GridOutFlags
{
+ DX::DX (const bool write_faces,
+ const bool write_all_faces) :
+ write_faces (write_faces),
+ write_all_faces (write_all_faces)
+ {};
+
+
+
Ucd::Ucd (const bool write_preamble,
const bool write_faces) :
write_preamble (write_preamble),
+void GridOut::set_flags (const GridOutFlags::DX &flags)
+{
+ dx_flags = flags;
+};
+
+
+
void GridOut::set_flags (const GridOutFlags::Ucd &flags)
{
ucd_flags = flags;
{
switch (output_format)
{
+ case dx:
+ return ".dx";
case gnuplot:
return ".gnuplot";
GridOut::OutputFormat
GridOut::parse_output_format (const std::string &format_name)
{
+ if (format_name == "dx")
+ return dx;
+
if (format_name == "ucd")
return ucd;
std::string GridOut::get_output_format_names ()
{
- return "gnuplot|eps|ucd";
+ return "dx|gnuplot|eps|ucd";
};
#include <cmath>
+template <int dim>
+void GridOut::write_dx (const Triangulation<dim> &tria,
+ std::ostream &out)
+{
+ AssertThrow (out, ExcIO());
+
+ // get the positions of the
+ // vertices and whether they are
+ // used.
+ const std::vector<Point<dim> > &vertices = tria.get_vertices();
+ const std::vector<bool> &vertex_used = tria.get_used_vertices();
+
+ const unsigned int n_vertices = tria.n_used_vertices();
+
+ typename Triangulation<dim>::active_cell_iterator cell=tria.begin_active();
+ const typename Triangulation<dim>::active_cell_iterator endc=tria.end();
+
+ // start with ucd data
+ out << n_vertices << ' '
+ << tria.n_active_cells() + (ucd_flags.write_faces ?
+ n_boundary_faces(tria) :
+ 0)
+ << " 0 0 0" // no data
+ << std::endl;
+
+ // actually write the vertices.
+ // note that we shall number them
+ // with first index 1 instead of 0
+ for (unsigned int i=0; i<vertices.size(); ++i)
+ if (vertex_used[i])
+ {
+ out << i+1 // vertex index
+ << " "
+ << vertices[i];
+ for (unsigned int d=dim+1; d<=3; ++d)
+ out << " 0"; // fill with zeroes
+ out << std::endl;
+ };
+
+ // write cells. Enumerate cells
+ // consecutively, starting with 1
+ unsigned int cell_index=1;
+ for (cell=tria.begin_active();
+ cell!=endc; ++cell, ++cell_index)
+ {
+ out << cell_index << ' '
+ << static_cast<unsigned int>(cell->material_id())
+ << " ";
+ switch (dim)
+ {
+ case 1: out << "line "; break;
+ case 2: out << "quad "; break;
+ case 3: out << "hex "; break;
+ default:
+ Assert (false, ExcNotImplemented());
+ };
+
+ // it follows a list of the
+ // vertices of each cell. in 1d
+ // this is simply a list of the
+ // two vertices, in 2d its counter
+ // clockwise, as usual in this
+ // library. in 3d, the same applies
+ // (special thanks to AVS for
+ // numbering their vertices in a
+ // way compatible to deal.II!)
+ //
+ // technical reference:
+ // AVS Developer's Guide, Release 4,
+ // May, 1992, p. E6
+ //
+ // note: vertex numbers are 1-base
+ for (unsigned int vertex=0; vertex<GeometryInfo<dim>::vertices_per_cell;
+ ++vertex)
+ out << cell->vertex_index(vertex)+1 << ' ';
+ out << std::endl;
+ };
+
+ // write faces with non-zero boundary
+ // indicator
+ if (ucd_flags.write_faces)
+ write_ucd_faces (tria, cell_index, out);
+
+ AssertThrow (out, ExcIO());
+};
+
+
+
template <int dim>
void GridOut::write_ucd (const Triangulation<dim> &tria,
std::ostream &out)
{
switch (output_format)
{
+ case dx:
+ write_dx (tria, out);
+ return;
+
case ucd:
write_ucd (tria, out);
return;