new grid output for OpenDX

git-svn-id: https://svn.dealii.org/trunk@5790 0785d39b-7218-0410-832d-ea1e28bc413d

diff --git a/deal.II/deal.II/include/grid/geometry_info.h b/deal.II/deal.II/include/grid/geometry_info.h
index b36a9c9b6b9324350f475bfbb2a16c5f1e55d14a..37bc3a79d9439177096f0f30912471fb2616d8e1 100644 (file)
--- a/deal.II/deal.II/include/grid/geometry_info.h
+++ b/deal.II/deal.II/include/grid/geometry_info.h
@@ -127,6 +127,28 @@ struct GeometryInfo<1>
                                      */
     static const unsigned int opposite_face[faces_per_cell];
     
+
+                                    /**
+                                     * Rearrange verices for OpenDX
+                                     * output.  For a cell being
+                                     * written in OpenDX format, each
+                                     * entry in this field contains
+                                     * the number of a vertex in
+                                     * @p{deal.II} that corresponds
+                                     * to the DX numbering at this
+                                     * location.
+                                     *
+                                     * Typical example: write a cell
+                                     * and arrange the vertices, such
+                                     * that OpenDX understands them.
+                                     *
+                                     * \begin{verbatim}
+                                     * for (i=0; i< n_vertices; ++i)
+                                     *   out << cell->vertex(dx_to_deal[i]);
+                                     * \end{verbatim}
+                                     */
+    static const unsigned int dx_to_deal[vertices_per_cell];
+
                                     /**
                                      * This field store which child cells
                                      * are adjacent to a certain face of
@@ -355,6 +377,27 @@ struct GeometryInfo<2>
                                      */
     static const unsigned int opposite_face[faces_per_cell];
     
+                                    /**
+                                     * Rearrange verices for OpenDX
+                                     * output.  For a cell being
+                                     * written in OpenDX format, each
+                                     * entry in this field contains
+                                     * the number of a vertex in
+                                     * @p{deal.II} that corresponds
+                                     * to the DX numbering at this
+                                     * location.
+                                     *
+                                     * Typical example: write a cell
+                                     * and arrange the vertices, such
+                                     * that OpenDX understands them.
+                                     *
+                                     * \begin{verbatim}
+                                     * for (i=0; i< n_vertices; ++i)
+                                     *   out << cell->vertex(dx_to_deal[i]);
+                                     * \end{verbatim}
+                                     */
+    static const unsigned int dx_to_deal[vertices_per_cell];
+
                                     /**
                                      * This field store which child cells
                                      * are adjacent to a certain face of
@@ -493,6 +536,27 @@ struct GeometryInfo<3>
                                      */
     static const unsigned int opposite_face[faces_per_cell];
     
+                                    /**
+                                     * Rearrange verices for OpenDX
+                                     * output.  For a cell being
+                                     * written in OpenDX format, each
+                                     * entry in this field contains
+                                     * the number of a vertex in
+                                     * @p{deal.II} that corresponds
+                                     * to the DX numbering at this
+                                     * location.
+                                     *
+                                     * Typical example: write a cell
+                                     * and arrange the vertices, such
+                                     * that OpenDX understands them.
+                                     *
+                                     * \begin{verbatim}
+                                     * for (i=0; i< n_vertices; ++i)
+                                     *   out << cell->vertex(dx_to_deal[i]);
+                                     * \end{verbatim}
+                                     */
+    static const unsigned int dx_to_deal[vertices_per_cell];
+
                                     /**
                                      * This field store which child cells
                                      * are adjacent to a certain face of

diff --git a/deal.II/deal.II/include/grid/grid_out.h b/deal.II/deal.II/include/grid/grid_out.h
index 3ef69a2a949cc119f0370245b150b23c93e56efc..b79105d73b239683e6eca50c0a023b524a7643ff 100644 (file)
--- a/deal.II/deal.II/include/grid/grid_out.h
+++ b/deal.II/deal.II/include/grid/grid_out.h
@@ -40,21 +40,38 @@ namespace GridOutFlags
   struct DX
   {
                                     /**
-                                     * Write faces instead of
-                                     * cells. Defaults to false.
+                                     * Write cells.
+                                     */
+    bool write_cells;
+      
+                                    /**
+                                     * Write faces.
                                      */
     bool write_faces;
+
+                                    /**
+                                     * Write field with diameters.
+                                     */
+    bool write_diameter;
+
+                                    /**
+                                     * Write field with area/volume.
+                                     */
+    bool write_measure;      
+                                     
                                     /**
                                      * Write all faces, including
-                                     * interior faces. Default is to
-                                     * write boundary faces only.
+                                     * interior faces.
                                      */
     bool write_all_faces;
     
                                     /**
                                      * Constructor.
                                      */
-    DX (const bool write_faces = false,
+    DX (const bool write_cells = true,
+       const bool write_faces = false,
+       const bool write_diameter = false,
+       const bool write_measure = false,
        const bool write_all_faces = false);  
   };                      
     
@@ -714,10 +731,10 @@ class GridOut
     DeclException0 (ExcInvalidState);
 
   private:
-                                  /**
-                                   * Flags for OpenDX output.
-                                   */
-  GridOutFlags::DX dx_flags;
+                                    /**
+                                     * Flags for OpenDX output.
+                                     */
+    GridOutFlags::DX dx_flags;
   
                                     /**
                                      * Flags for UCD output. Can be

diff --git a/deal.II/deal.II/source/grid/geometry_info.cc b/deal.II/deal.II/source/grid/geometry_info.cc
index b7e756c8c3d5083bd94836ec43d878ff8b3da794..812106d8d72778d029e88aadd3167bfe98babc13 100644 (file)
--- a/deal.II/deal.II/source/grid/geometry_info.cc
+++ b/deal.II/deal.II/source/grid/geometry_info.cc
@@ -66,6 +66,19 @@ const unsigned int GeometryInfo<3>::opposite_face[GeometryInfo<3>::faces_per_cel
 
 
 
+const unsigned int GeometryInfo<1>::dx_to_deal[GeometryInfo<1>::vertices_per_cell]
+= { 0, 1};
+
+
+const unsigned int GeometryInfo<2>::dx_to_deal[GeometryInfo<2>::vertices_per_cell]
+= { 0, 3, 1, 2};
+
+
+const unsigned int GeometryInfo<3>::dx_to_deal[GeometryInfo<3>::vertices_per_cell]
+= { 0, 3, 4, 7, 1, 2, 5, 6};
+
+
+
 
 
 unsigned int

diff --git a/deal.II/deal.II/source/grid/grid_out.all_dimensions.cc b/deal.II/deal.II/source/grid/grid_out.all_dimensions.cc
index 40e605e5f1a494f37365f68a7e398ea12577df93..b98eab480cdb28bc9ac6ec89e1a79bbad07545aa 100644 (file)
--- a/deal.II/deal.II/source/grid/grid_out.all_dimensions.cc
+++ b/deal.II/deal.II/source/grid/grid_out.all_dimensions.cc
@@ -17,9 +17,15 @@
 
 namespace GridOutFlags
 {
-  DX::DX (const bool write_faces,
+  DX::DX (const bool write_cells,
+         const bool write_faces,
+         const bool write_diameter,
+         const bool write_measure,
          const bool write_all_faces) :
+    write_cells (write_cells),
     write_faces (write_faces),
+    write_diameter (write_diameter),
+    write_measure (write_measure),
     write_all_faces (write_all_faces)
   {};
 

diff --git a/deal.II/deal.II/source/grid/grid_out.cc b/deal.II/deal.II/source/grid/grid_out.cc
index d353603b12607d5f518e496a8e4ad2a07badfb03..6dec185165ed7f736b1484a008038d9d6783ee1e 100644 (file)
--- a/deal.II/deal.II/source/grid/grid_out.cc
+++ b/deal.II/deal.II/source/grid/grid_out.cc
@@ -27,14 +27,232 @@
 #include <cmath>
 
 
+#if deal_II_dimension == 1
+
+
+template <int dim>
+void GridOut::write_dx (const Triangulation<dim> &,
+                       std::ostream             &)
+{
+  Assert (false, ExcNotImplemented());
+}
+
+#else
+
+
 template <int dim>
-void GridOut::write_dx (const Triangulation<dim> &/*tria*/,
+void GridOut::write_dx (const Triangulation<dim> &tria,
                        std::ostream             &out) 
 {
-  Assert(false, ExcNotImplemented());
+//  Assert(false, ExcNotImplemented());
   AssertThrow (out, ExcIO());
-};
+                                  // Copied and adapted from write_ucd
+  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;
+  const typename Triangulation<dim>::active_cell_iterator endc=tria.end();
+
+  
+                                  // write the vertices
+  out << "object \"vertices\" class array type float rank 1 shape " << dim
+      << " items " << n_vertices << " data follows"
+      << std::endl;
+  
+  for (unsigned int i=0; i<vertices.size(); ++i)
+    if (vertex_used[i])
+      {
+       out << '\t' << vertices[i] << std::endl;
+      };
+  
+                                  // write cells or faces
+  const bool write_cells = dx_flags.write_cells;
+  const bool write_faces = (dim>1) ? dx_flags.write_faces : false;
+  
+  const unsigned int n_cells = tria.n_active_cells();
+  const unsigned int n_faces = tria.n_active_cells()
+                              * GeometryInfo<dim>::faces_per_cell;
+
+  const unsigned int n_vertices_per_cell = GeometryInfo<dim>::vertices_per_cell;
+  const unsigned int n_vertices_per_face = GeometryInfo<dim>::vertices_per_face;
+  
+  if (write_cells)
+    {
+      out << "object \"cells\" class array type int rank 1 shape "
+         << n_vertices_per_cell
+         << " items " << n_cells << " data follows" << std::endl;
+      
+      for (cell = tria.begin_active(); cell != endc; ++cell)
+       {
+         for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
+           out << '\t' << cell->vertex_index(GeometryInfo<dim>::dx_to_deal[v]);
+         out << std::endl;
+       }
+      out << "attribute \"element type\" string \"";
+      if (dim==1) out << "lines";
+      if (dim==2) out << "quads";
+      if (dim==3) out << "cubes";
+      out << "\"" << std::endl
+         << "attribute \"ref\" string \"positions\"" << std::endl << std::endl;
+
+                                      // Additional cell information
+      
+      out << "object \"material\" class array type int rank 0 items "
+         << n_cells << " data follows" << std::endl;
+      for (cell = tria.begin_active(); cell != endc; ++cell)
+       out << ' ' << (unsigned int)cell->material_id();
+      out  << std::endl
+          << "attribute \"dep\" string \"connections\"" << std::endl << std::endl;
+      
+      out << "object \"level\" class array type int rank 0 items "
+         << n_cells << " data follows" << std::endl;
+      for (cell = tria.begin_active(); cell != endc; ++cell)
+       out << ' ' << cell->level();
+      out  << std::endl
+          << "attribute \"dep\" string \"connections\"" << std::endl << std::endl;
+      
+      if (dx_flags.write_measure)
+       {
+         out << "object \"measure\" class array type float rank 0 items "
+             << n_cells << " data follows" << std::endl;
+         for (cell = tria.begin_active(); cell != endc; ++cell)
+           out << '\t' << cell->measure();
+         out  << std::endl
+              << "attribute \"dep\" string \"connections\"" << std::endl << std::endl;
+       }
+      
+      if (dx_flags.write_diameter)
+       {
+         out << "object \"diameter\" class array type float rank 0 items "
+             << n_cells << " data follows" << std::endl;
+         for (cell = tria.begin_active(); cell != endc; ++cell)
+           out << '\t' << cell->diameter();
+         out  << std::endl
+              << "attribute \"dep\" string \"connections\"" << std::endl << std::endl;
+       }
+    }
+  
+  if (write_faces)
+    {
+      out << "object \"faces\" class array type int rank 1 shape "
+         << n_vertices_per_face
+         << " items " << n_faces << " data follows"
+         << std::endl;
+
+      for (cell = tria.begin_active(); cell != endc; ++cell)
+       {
+         for (unsigned int f=0;f<GeometryInfo<dim>::faces_per_cell;++f)
+           {
+             typename Triangulation<dim>::face_iterator face = cell->face(f);
+             
+             for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_face; ++v)
+               out << '\t' << face->vertex_index(GeometryInfo<dim-1>::dx_to_deal[v]);
+             out << std::endl;
+           }
+       }
+      out << "attribute \"element type\" string \"";
+      if (dim==2) out << "lines";
+      if (dim==3) out << "quads";
+      out << "\"" << std::endl
+         << "attribute \"ref\" string \"positions\"" << std::endl << std::endl;
+      
+
+                                      // Additional face information
+      
+      out << "object \"boundary\" class array type int rank 0 items "
+         << n_faces << " data follows" << std::endl;
+      for (cell = tria.begin_active(); cell != endc; ++cell)
+       {
+                                          // Little trick to get -1
+                                          // for the interior
+         for (unsigned int f=0;f<GeometryInfo<dim>::faces_per_cell;++f)
+           out << ' ' << (int)(signed char)cell->face(f)->boundary_indicator();
+         out << std::endl;
+       }
+      out << "attribute \"dep\" string \"connections\"" << std::endl << std::endl;
+      
+      if (dx_flags.write_measure)
+       {
+         out << "object \"face measure\" class array type float rank 0 items "
+             << n_faces << " data follows" << std::endl;
+         for (cell = tria.begin_active(); cell != endc; ++cell)
+           {
+             for (unsigned int f=0;f<GeometryInfo<dim>::faces_per_cell;++f)
+               out << ' ' << cell->face(f)->measure();
+             out << std::endl;
+           }
+         out << "attribute \"dep\" string \"connections\"" << std::endl << std::endl;
+       }
+
+      if (dx_flags.write_diameter)
+       {
+         out << "object \"face diameter\" class array type float rank 0 items "
+             << n_faces << " data follows" << std::endl;
+         for (cell = tria.begin_active(); cell != endc; ++cell)
+           {
+             for (unsigned int f=0;f<GeometryInfo<dim>::faces_per_cell;++f)
+               out << ' ' << cell->face(f)->diameter();
+             out << std::endl;
+           }
+         out << "attribute \"dep\" string \"connections\"" << std::endl << std::endl;
+       }
+
+    }
+  
+
+                                  // Write additional face information
+  
+  if (write_faces)
+    {
+      
+    }
+  else
+    {
+     }
+
+                                  // The wrapper
+  out << "object \"deal data\" class field" << std::endl
+      << "component \"positions\" value \"vertices\"" << std::endl
+      << "component \"connections\" value \"cells\"" << std::endl;
+
+  if (write_cells)
+    {
+      out << "object \"cell data\" class field" << std::endl
+         << "component \"positions\" value \"vertices\"" << std::endl
+         << "component \"connections\" value \"cells\"" << std::endl;
+      out << "component \"material\" value \"material\"" << std::endl;
+      out << "component \"level\" value \"level\"" << std::endl;
+      if (dx_flags.write_measure)
+       out << "component \"measure\" value \"measure\"" << std::endl;
+      if (dx_flags.write_diameter)
+       out << "component \"diameter\" value \"diameter\"" << std::endl;
+    }
 
+  if (write_faces)
+    {
+      out << "object \"face data\" class field" << std::endl
+         << "component \"positions\" value \"vertices\"" << std::endl
+         << "component \"connections\" value \"faces\"" << std::endl;
+      out << "component \"boundary\" value \"boundary\"" << std::endl;
+      if (dx_flags.write_measure)
+       out << "component \"measure\" value \"face measure\"" << std::endl;
+      if (dx_flags.write_diameter)
+       out << "component \"diameter\" value \"face diameter\"" << std::endl;
+    }
+  
+  out << std::endl
+      << "object \"grid data\" class group" << std::endl;
+    if (write_cells)
+      out << "member \"cells\" value \"cell data\"" << std::endl;
+    if (write_faces)
+      out << "member \"faces\" value \"face data\"" << std::endl;
+  out << "end" << std::endl;  
+}
+
+
+#endif
 
 
 template <int dim>