]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
dx output declared
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 30 Nov 2001 07:33:59 +0000 (07:33 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 30 Nov 2001 07:33:59 +0000 (07:33 +0000)
git-svn-id: https://svn.dealii.org/trunk@5316 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/grid/grid_out.h
deal.II/deal.II/source/grid/grid_out.all_dimensions.cc
deal.II/deal.II/source/grid/grid_out.cc

index e98fab327e25d319d627d859df345ab73273490c..a000405a5ca84dd62a35d40103c2b1a53d08d8e9 100644 (file)
@@ -34,6 +34,30 @@ template <int dim> class Mapping;
  */
 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.
@@ -386,7 +410,19 @@ class GridOut
                                      * 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
@@ -583,34 +619,37 @@ class GridOut
                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
@@ -618,6 +657,7 @@ class GridOut
                                      * 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}.
@@ -674,10 +714,14 @@ class GridOut
     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;
index 3ac37f64c9bde209566cadbae63b73ed7898b32e..cc01ce2146d3224da1868912241b160bc8a3c68f 100644 (file)
 
 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),
@@ -91,6 +99,13 @@ namespace GridOutFlags
 
 
 
+void GridOut::set_flags (const GridOutFlags::DX &flags) 
+{
+  dx_flags = flags;
+};
+
+
+
 void GridOut::set_flags (const GridOutFlags::Ucd &flags) 
 {
   ucd_flags = flags;
@@ -131,6 +146,8 @@ GridOut::default_suffix (const OutputFormat output_format)
 {
   switch (output_format) 
     {
+      case dx:
+        return ".dx";
       case gnuplot: 
            return ".gnuplot";
 
@@ -151,6 +168,9 @@ GridOut::default_suffix (const OutputFormat output_format)
 GridOut::OutputFormat
 GridOut::parse_output_format (const std::string &format_name)
 {
+  if (format_name == "dx")
+    return dx;
+
   if (format_name == "ucd")
     return ucd;
 
@@ -169,7 +189,7 @@ GridOut::parse_output_format (const std::string &format_name)
 
 std::string GridOut::get_output_format_names () 
 {
-  return "gnuplot|eps|ucd";
+  return "dx|gnuplot|eps|ucd";
 };
 
 
index 4329eb492cc8e14c993aeb81e8018028464e00a4..968a2ffdd46bd73dc08e45b0b28221c82d44c7e4 100644 (file)
 #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) 
@@ -806,6 +894,10 @@ void GridOut::write (const Triangulation<dim> &tria,
 {
   switch (output_format)
     {
+      case dx:
+           write_dx (tria, out);
+           return;
+
       case ucd:
            write_ucd (tria, out);
            return;

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.