]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implement a way to pass flags to the output functions for GridOut.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 29 Jun 1999 18:47:27 +0000 (18:47 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 29 Jun 1999 18:47:27 +0000 (18:47 +0000)
git-svn-id: https://svn.dealii.org/trunk@1507 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 8e66cf8a05fa6731a3ed8dbf13ac047b96ab7145..05b63676c2a970248373e7e66bc84132a902bfff 100644 (file)
@@ -13,7 +13,8 @@
 /**
  * This class provides a means to output a triangulation to a file in different
  * formats. Presently provided are functions to write it in GNUPLOT and
- * encapsulated postscript format. There are also function to dispatch to
+ * encapsulated postscript format, as well as in UCD format.
+ * There are also function to dispatch to
  * the different output functions based on a parameter given, e.g., through
  * a parameter file, thus making user programs invariant of the number and
  * names of the file formats implemented in this class. The main advantage of
  * write the pure geometrical information involved here.
  *
  *
+ * \subsection{UCD format}
+ * UCD (unstructured cell data) is the format used by AVS and some other
+ * programs. It is described in the AVS developer's guide. Besides the usual
+ * output of the grid only, you can decide through additional flags (see
+ * below, and the documentation of the #GridOut::UcdFlags# class) whether
+ * boundary faces with non-zero boundary indicator shall be written to the
+ * file explicitely. This is useful, if you want to re-read the grid later
+ * on, since #deal.II# sets the boundary indicator to zero by default;
+ * therefore, to obtain the same triangulation as before, you have to specify
+ * faces with differing boundary indicators explicitely, which is done by
+ * this flag.
+ *
+ * Names and values of further flags controlling the output can be found
+ * in the documentation of the #GridOut::UcdFlags# class.
+ *
+ *
  * \subsection{GNUPLOT format}
  * In GNUPLOT format, each cell is written as a sequence of lines. There is not
  * much more to be said here since the actual representation as well as the
  * viewing angle (for 3d plots) is done by GNUPLOT itself.
  *
- * One additional feature, however, is worth being mentioned: after each point
- * denoting one of the end points of a line, the level of the respective cell
- * is written. Therefore, if you let GNUPLOT draw a 2d grid as a 3d plot, you
- * will see more refined cells being raised against cells with less refinement.
- * Also, if you draw a cut through a 3d grid, you can extrude the refinement
- * level in the direction orthogonal to the cut plane.
+ * One additional feature, however, is worth being mentioned: if
+ * switched on by using the flags controlling GNUPLOT output (see
+ * below on how to do this), after each point denoting one of the end
+ * points of a line, the level of the respective cell is
+ * written. Therefore, if you let GNUPLOT draw a 2d grid as a 3d plot,
+ * you will see more refined cells being raised against cells with
+ * less refinement.  Also, if you draw a cut through a 3d grid, you
+ * can extrude the refinement level in the direction orthogonal to the
+ * cut plane.
  *
  * A more useful application of this feature is the following: if you use the
  * GNUPLOT command (for a 2d grid here)
  * plot of those cells only that belong to level 3 in this example. This
  * way, it is easy to produce plots of the different levels of grid.
  *
+ * Names and values of additional flags controlling the output can be found
+ * in the documentation of the #GridOut::GnuplotFlags# class.
+ *
  *
  * \subsection{Encapsulated postscript format}
  * In this format, each line of the triangulation is written separately. We
- * scale the picture such that x-values range between zero and 300 and the y-range
- * is scaled with the same factor. The bounding box is close to the triangulation
- * on all four sides, without an extra frame. The line width is chosen to be 0.5,
- * which is a thin line relative to the extension of the picture of 300.
+ * scale the picture such that either x-values or y-values range between zero and
+ * a fixed size. The other axis is scaled by the same factor. Which axis is
+ * taken to compute the scale and the size of the box it shall fit into is
+ * determined by the output flags (see below, and the documentation of the
+ * #GridOut::EpsFlags# class).
+ *
+ * The bounding box is close to the triangulation on all four sides, without an
+ * extra frame. The line width is chosen to be 0.5 by default, but can be
+ * changed. The line width is to be compared with the extension of the picture,
+ * of which the default is 300.
+ *
+ * Names and values of additional flags controlling the output can be found
+ * in the documentation of the #GridOut::GnuplotFlags# class. Especially
+ * the viewpoint for three dimensional grids is of importance here.
  * 
  *
  * \subsection{Usage}
  * Usage is simple: either you use the direct form
  * \begin{verbatim}
  *   ofstream output_file("some_filename");
- *   GridOut::write_gnuplot (tria, output_file);
+ *   GridOut().write_gnuplot (tria, output_file);
  * \end{verbatim}
  * if you know which format you want to have, or if you want the format to be
  * a runtime parameter, you can write
  *   GridOut::OutputFormat grid_format =
  *                   GridOut::parse_output_format(get_format_name_from_somewhere());
  *   ofstream output_file("some_filename" + GridOut::default_suffix(output_format));
- *   GridOut::write (tria, output_file, output_format);
+ *   GridOut().write (tria, output_file, output_format);
  * \end{verbatim}
  * The function #get_output_format_names()# provides a list of possible names of
  * output formats in a string that is understandable by the #ParameterHandler# class.
  *
+ * Note that here, we have created an unnamed object of type #GridOut# and called
+ * one of its #write_*# functions. This looks like as if the respective function
+ * could really be made #static#. This was not done in order to allow for
+ * parameters to be passed to the different output functions in a way compatible
+ * with the scheme of allowing the right output format to be selected at run-time
+ * through the generic #write# function.
+ *
+ * In order to explain this, consider each function had one or more additional
+ * parameters giving the details of output, for example position of the spectator
+ * for 3d meshed, line thicknesses, etc. While this would allow each output
+ * function any flexibility it needs, it would not allow us to use the generic
+ * function #write# which is given a parameter determining the output format,
+ * since it is impractical to give it a list of parameters for each and every
+ * output format supported which it may then pass on to the respective output
+ * function.
+ *
+ * Rather, we have chosen to let each object of this class #GridOut# have a
+ * set of parameters for each supported output format. These are collected
+ * in structures #EpsFlags#, #GnuplotFlags#, etc and you can set your preferred
+ * flags like this:
+ * \begin{verbatim}
+ *   GridOut grid_out;
+ *   GridOut::UcdFlags ucd_flags;
+ *   ...    // set some fields in ucd_flags
+ *   grid_out.set_flags (ucd_flags);
+ *   ...
+ *   ...    // write some file with data_out
+ * \end{verbatim}
+ * The respective output function then use the so-set flags. By default, they
+ * are set to reasonable values as described above and in the documentation
+ * of the different flags structures. Resetting the flags can
+ * be done by calling #grid_out.set_flags (GridOut::UcdFlags());#, since the
+ * default constructor of each of the flags structures sets the parameters
+ * to their initial values.
  *
- * @author Wolfgang Bangerth, 1999; postscript format based on an implementation by Stefan Nauber, 1999
- */
+ * The advantage of this approach is that it is possible to change the flags
+ * of one or more output formats according to your needs and later use the
+ * generic #write# function; the actual output function then called will
+ * use the flags as set before.
+ *
+ * Note that some of the structures describing the flags of the different
+ * output formats are empty since the respective format does not support
+ * any flags. The structure and the #set_flags# function are provided
+ * anyway. Note also that some of the structures may differ between the
+ * dimensions supported by this class; they then have a template parameter,
+ * as usual.
+ *
+ *
+ * @author Wolfgang Bangerth, 1999; postscript format based on an implementation by Stefan Nauber, 1999 */
 class GridOut 
 {
   public:
+
+                                    /**
+                                     * Flags describing the details of
+                                     * output in UCD format.
+                                     */
+    struct UcdFlags 
+    {
+                                        /**
+                                         * Write a comment at the beginning
+                                         * of the file stating the date of
+                                         * creation and some other data.
+                                         * While this is supported by the
+                                         * UCD format (and the AVS program),
+                                         * some other programs get confused
+                                         * by this, so you can switch it off
+                                         * this way.
+                                         *
+                                         * Default: #true#.
+                                         */
+       bool write_preamble;
+       
+                                        /**
+                                         * When writing a mesh, write boundary
+                                         * faces explicitely if their boundary
+                                         * indicator is not the default
+                                         * boundary indicator, which is zero.
+                                         * This is necessary if you later
+                                         * want to re-read the grid and want
+                                         * to get the same boundary indicators
+                                         * for the different parts of the
+                                         * boundary of the triangulation.
+                                         *
+                                         * It is not necessary if you only want
+                                         * to write the triangulation to
+                                         * view or print it.
+                                         *
+                                         * Default: #false#.
+                                         */
+       bool write_faces;
+
+                                        /**
+                                         * Constructor.
+                                         */
+       UcdFlags (const bool write_preamble = true,
+                 const bool write_faces    = false);
+    };
+
+
+                                    /**
+                                     * Flags describing the details of
+                                     * output in GNUPLOT format.
+                                     */
+    struct GnuplotFlags 
+    {
+                                        /**
+                                         * Write the number of each cell into
+                                         * the output file before starting
+                                         * with the lines it is composed of,
+                                         * as a comment. This might be useful
+                                         * if you want to find out details
+                                         * about the grid, for example the
+                                         * position of cells of which you
+                                         * know the number. It enlarges
+                                         * the size of the output
+                                         * significantly, however.
+                                         *
+                                         * Default: #false#.
+                                         */
+       bool write_cell_numbers;
+
+                                        /**
+                                         * Write the level of a cell as
+                                         * an additional column after the
+                                         * coordinates. See the general
+                                         * documentation of this class
+                                         * for an example of the use
+                                         * of this feature.
+                                         *
+                                         * Default: #false#.
+                                         */
+       bool write_level;
+       
+                                        /**
+                                         * Constructor.
+                                         */
+       GnuplotFlags (const bool write_cell_number = false,
+                     const bool write_level       = false);
+    };
+
+                                    /**
+                                     * Flags describing the details of
+                                     * output for encapsulated postscript.
+                                     * In this structure, the flags common
+                                     * to all dimensions are listed. Flags
+                                     * which are specific to one space
+                                     * dimension only are listed in derived
+                                     * classes.
+                                     *
+                                     * By default, the size of the picture
+                                     * is scaled such that the width equals
+                                     * 300 units.
+                                     */
+    struct EpsFlagsBase
+    {
+                                        /**
+                                         * Enum denoting the possibilities
+                                         * whether the scaling should be done
+                                         * such that the given #size# equals
+                                         * the width or the height of
+                                         * the resulting picture.
+                                         */
+       enum SizeType {
+             width, height
+       };
+
+                                        /**
+                                         * See above. Default is #width#.
+                                         */
+       SizeType size_type;
+       
+                                        /**
+                                         * Width or height of the output
+                                         * as given in postscript units
+                                         * This usually is given by the
+                                         * strange unit 1/72 inch. Whether
+                                         * this is height or width is
+                                         * specified by the flag
+                                         * #size_type#.
+                                         *
+                                         * Default is 300.
+                                         */
+       unsigned int size;
+
+                                        /**
+                                         * Width of a line in postscript
+                                         * units. Default is 0.5.
+                                         */
+       double line_width;
+       
+                                        /**
+                                         * Constructor.
+                                         */
+       EpsFlagsBase (const SizeType     size_type  = width,
+                     const unsigned int size       = 300,
+                     const double       line_width = 0.5);
+    };
+    
+       
+                                    /**
+                                     * Flags describing the details of
+                                     * output for encapsulated postscript
+                                     * for all dimensions not explicitely
+                                     * specialized below. Some flags that
+                                     * are common to all dimensions are
+                                     * listed in the base class
+                                     */
+    template <int dim>
+    struct EpsFlags : public EpsFlagsBase 
+    {};
+
+                                    /**
+                                     * Flags specific to the output of
+                                     * grids in three space dimensions.
+                                     */
+    template <>
+    struct EpsFlags<3> : public EpsFlagsBase 
+    {
+                                        /**
+                                         * Angle of the line origin-viewer
+                                         * against the z-axis in degrees
+                                         *
+                                         * Default is the Gnuplot-default
+                                         * of 60.
+                                         */
+       double azimut_angle;
+
+                                        /**
+                                         * Angle by which the viewers
+                                         * position projected onto the
+                                         * x-y-plane is rotated around
+                                         * the z-axis, in positive sense
+                                         * when viewed from above. The
+                                         * unit are degrees, and zero
+                                         * equals a position above or below
+                                         * the negative y-axis.
+                                         *
+                                         * Default is the Gnuplot-default
+                                         * of 30.
+                                         */
+       double turn_angle;
+
+                                        /**
+                                         * Constructor.
+                                         */
+       EpsFlags (const double        azimut_angle    = 60,
+                 const double        turn_angle      = 30) :
+                       azimut_angle (azimut_angle),
+                       turn_angle (turn_angle)      {};
+    };
+    
+       
+    
                                     /**
                                      * Declaration of a name for each of the
                                      * different output formats.
@@ -123,7 +404,37 @@ class GridOut
     void write (const Triangulation<dim> &tria,
                ostream                  &out,
                const OutputFormat        output_format);
-    
+
+                                    /**
+                                     * Set the flags to be used for output
+                                     * in UCD format.
+                                     */
+    void set_flags (const UcdFlags &ucd_flags);
+
+                                    /**
+                                     * Set the flags to be used for output
+                                     * in GNUPLOT format.
+                                     */
+    void set_flags (const GnuplotFlags &gnuplot_flags);
+
+                                    /**
+                                     * Set the flags to be used for output
+                                     * in 1d EPS output.
+                                     */
+    void set_flags (const EpsFlags<1> &eps_flags);
+
+                                    /**
+                                     * Set the flags to be used for output
+                                     * in 2d EPS output.
+                                     */
+    void set_flags (const EpsFlags<2> &eps_flags);
+
+                                    /**
+                                     * Set the flags to be used for output
+                                     * in 3d EPS output.
+                                     */
+    void set_flags (const EpsFlags<3> &eps_flags);
+
                                     /**
                                      * Provide a function which tells us which
                                      * suffix with a given output format
@@ -185,6 +496,46 @@ class GridOut
 
   private:
 
+                                    /**
+                                     * Flags to be used upon output of UCD
+                                     * data. Can be changed by using the
+                                     * #set_flags# function.
+                                     */
+    UcdFlags     ucd_flags;
+
+                                    /**
+                                     * Flags to be used upon output of GNUPLOT
+                                     * data. Can be changed by using the
+                                     * #set_flags# function.
+                                     */
+    GnuplotFlags gnuplot_flags;
+
+                                    /**
+                                     * Flags to be used upon output of EPS
+                                     * data in one space dimension. Can be
+                                     * changed by using the #set_flags#
+                                     * function.
+                                     */
+    EpsFlags<1>  eps_flags_1;
+
+                                    /**
+                                     * Flags to be used upon output of EPS
+                                     * data in two space dimensions. Can be
+                                     * changed by using the #set_flags#
+                                     * function.
+                                     */
+    EpsFlags<2>  eps_flags_2;
+
+                                    /**
+                                     * Flags to be used upon output of EPS
+                                     * data in three space dimensions. Can be
+                                     * changed by using the #set_flags#
+                                     * function.
+                                     */
+    EpsFlags<3>  eps_flags_3;
+    
+    
+
                                     /**
                                      * Write the grid information about
                                      * faces to #out#. Only those faces
@@ -238,6 +589,7 @@ class GridOut
 
 
 
+
 /*----------------------------   grid_out.h     ---------------------------*/
 /* end of #ifndef __grid_out_H */
 #endif
index 99610be3839e8be43b50fd6b99dc2caab191a4cc..77651cc8f3d169aa608ff814d7e96312297cf436 100644 (file)
@@ -3,6 +3,79 @@
 #include <basic/grid_out.h>
 
 
+GridOut::UcdFlags::UcdFlags (const bool write_preamble,
+                            const bool write_faces) :
+               write_preamble (write_preamble),
+               write_faces (write_faces)
+{};
+
+
+
+GridOut::GnuplotFlags::GnuplotFlags (const bool write_cell_numbers,
+                                    const bool write_level) :
+               write_cell_numbers (write_cell_numbers),
+               write_level (write_level)
+{};
+
+
+
+GridOut::EpsFlagsBase::EpsFlagsBase (const SizeType     size_type,
+                                    const unsigned int size,
+                                    const double       line_width) :
+               size_type (size_type),
+               size (size),
+               line_width (line_width)
+{};
+
+
+
+// egcs 1.1.2 does not understand this, so I made it inlined in the
+// class declaration
+//
+// template <>
+// GridOut::EpsFlags<3>::EpsFlags (const double        azimut_angle,
+//                             const double        turn_angle) :
+//             azimut_angle (azimut_angle),
+//             turn_angle (turn_angle)
+// {};
+
+
+
+void GridOut::set_flags (const UcdFlags &flags) 
+{
+  ucd_flags = flags;
+};
+
+
+
+void GridOut::set_flags (const GnuplotFlags &flags) 
+{
+  gnuplot_flags = flags;
+};
+
+
+
+void GridOut::set_flags (const EpsFlags<1> &flags) 
+{
+  eps_flags_1 = flags;
+};
+
+
+
+void GridOut::set_flags (const EpsFlags<2> &flags) 
+{
+  eps_flags_2 = flags;
+};
+
+
+
+void GridOut::set_flags (const EpsFlags<3> &flags) 
+{
+  eps_flags_3 = flags;
+};
+
+
+
 string GridOut::default_suffix (const OutputFormat output_format) 
 {
   switch (output_format) 
index 9d3ce9bb21b7b4a8782dbc5e994cb8a71914c860..de39f71c68f99dbe1b8133e9ca9fcc6f2e62355e 100644 (file)
@@ -45,9 +45,8 @@ void GridOut::write_ucd (const Triangulation<dim> &tria,
   
 
                                   // write preamble
-  if (true)
+  if (ucd_flags.write_preamble)
     {
-/*      
                                       // block this to have local
                                       // variables destroyed after
                                       // use
@@ -66,12 +65,13 @@ void GridOut::write_ucd (const Triangulation<dim> &tria,
          << "# For a description of the UCD format see the AVS Developer's guide."
          << endl
          << "#" << endl;
-*/
     };
 
                                   // start with ucd data
   out << n_vertices << ' '
-      << tria.n_active_cells() + n_boundary_faces(tria)
+      << tria.n_active_cells() + (ucd_flags.write_faces ?
+                                 n_boundary_faces(tria) :
+                                 0)
       << " 0 0 0"                  // no data
       << endl;
 
@@ -130,7 +130,8 @@ void GridOut::write_ucd (const Triangulation<dim> &tria,
 
                                   // write faces with non-zero boundary
                                   // indicator
-  write_ucd_faces (tria, cell_index, out);
+  if (ucd_flags.write_faces)
+    write_ucd_faces (tria, cell_index, out);
     
   AssertThrow (out, ExcIO());
 };
@@ -223,7 +224,8 @@ void GridOut::write_gnuplot (const Triangulation<dim> &tria,
   const typename Triangulation<dim>::active_cell_iterator  endc=tria.end();
   for (; cell!=endc; ++cell)
     {
-      out << "# cell " << cell << endl;
+      if (gnuplot_flags.write_cell_numbers)
+       out << "# cell " << cell << endl;
       
       switch (dim)
        {
@@ -286,7 +288,19 @@ void GridOut::write_eps (const Triangulation<dim> &tria,
                         ostream                  &out) 
 {
   typedef list<pair<Point<2>,Point<2> > > LineList;
-  
+
+                                  // get a pointer to the flags
+                                  // common to all dimensions,
+                                  // in order to avoid the recurring
+                                  // distinctions between
+                                  // eps_flags_1, eps_flags_2, ...
+  const EpsFlagsBase &eps_flags_base = (dim==1 ?
+                                       eps_flags_1 :
+                                       (dim==2 ?
+                                        eps_flags_2 :
+                                        (dim==3 ?
+                                         eps_flags_3 :
+                                         *(EpsFlagsBase*)0)));
   
   AssertThrow (out, ExcIO());
 
@@ -345,9 +359,15 @@ void GridOut::write_eps (const Triangulation<dim> &tria,
                                         // spectator to the origin.
                                         //
                                         // we chose here the viewpoint as in
-                                        // gnuplot
-       const double z_angle    = 60;
-       const double turn_angle = 30;
+                                        // gnuplot as default.
+                                        //
+                                        // note: the following might be wrong
+                                        // if one of the base vectors below
+                                        // is in direction of the viewer, but
+                                        // I am too tired at present to fix
+                                        // this
+       const double z_angle    = eps_flags_3.azimut_angle;
+       const double turn_angle = eps_flags_3.turn_angle;
        const double pi = 3.1415926536;
        const Point<dim> view_direction(-sin(z_angle * 2.*pi / 360.) * sin(turn_angle * 2.*pi / 360.),
                                        +sin(z_angle * 2.*pi / 360.) * cos(turn_angle * 2.*pi / 360.),
@@ -415,7 +435,10 @@ void GridOut::write_eps (const Triangulation<dim> &tria,
                                   // don't scale in y-direction to
                                   // preserve the shape of the
                                   // triangulation
-  const double scale = 300. / (x_max - x_min);
+  const double scale = (eps_flags_base.size /
+                       (eps_flags_base.size_type==EpsFlagsBase::width ?
+                        x_max - x_min :
+                        y_min - y_max));
   
   
                                   // now write preamble
@@ -440,7 +463,9 @@ void GridOut::write_eps (const Triangulation<dim> &tria,
                                         // lower left corner
          << "0 0 "
                                         // upper right corner
-         << "300 " << static_cast<unsigned int>( (y_max-y_min) * scale )
+         << static_cast<unsigned int>( (x_max-x_min) * scale )
+         << ' '
+         << static_cast<unsigned int>( (y_max-y_min) * scale )
          << endl;
 
                                       // define some abbreviations to keep
@@ -454,7 +479,7 @@ void GridOut::write_eps (const Triangulation<dim> &tria,
          << endl;
 
                                       // set fine lines
-      out << "0.5 setlinewidth" << endl;
+      out << eps_flags_base.line_width << " setlinewidth" << endl;
     };
 
                                   // now write the lines

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.