]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add Vtk to list of supported file formats.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 19 Mar 2001 10:03:41 +0000 (10:03 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 19 Mar 2001 10:03:41 +0000 (10:03 +0000)
git-svn-id: https://svn.dealii.org/trunk@4242 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/include/base/data_out_base.h
deal.II/base/source/data_out_base.all_dimensions.cc
deal.II/base/source/data_out_base.cc
deal.II/doc/news/2001/c-3-1.html
deal.II/doc/readme.html

index 05465b3fa0c2488970450cfe9a1c3d36344382cf..0c102d721c0e715055c7503d9bf799e40ec1882f 100644 (file)
@@ -266,6 +266,14 @@ class ParameterHandler;
  * (bi-,tri-)linear elements.
  *
  *
+ * @sect3{VTK format}
+ *
+ * This is the file format used by the Visualization Toolkit (VTK, see
+ * @url{http://www.kitware.com/vtk.html}), as
+ * described in their manual, section 14.3. It is similar to the GMV
+ * format, see there for more information.
+ *
+ *
  * @sect2{Output parameters}
  *
  * All functions take a parameter which is a structure of type @p{XFlags}, where
@@ -945,6 +953,63 @@ class DataOutBase
                                      * present no flags are implemented.
                                      */
     struct GmvFlags 
+    {
+      private:
+                                        /**
+                                         * Dummy entry to suppress compiler
+                                         * warnings when copying an empty
+                                         * structure. Remove this member
+                                         * when adding the first flag to
+                                         * this structure (and remove the
+                                         * @p{private} as well).
+                                         */
+       int dummy;
+
+      public:
+                                        /**
+                                         * Declare all flags with name
+                                         * and type as offered by this
+                                         * class, for use in input files.
+                                         */
+       static void declare_parameters (ParameterHandler &prm);
+
+                                        /**
+                                         * Read the parameters declared in
+                                         * @p{declare_parameters} and set the
+                                         * flags for this output format
+                                         * accordingly.
+                                         *
+                                         * The flags thus obtained overwrite
+                                         * all previous contents of this object.
+                                         */
+       void parse_parameters (ParameterHandler &prm);
+
+                                        /**
+                                         * Determine an estimate for
+                                         * the memory consumption (in
+                                         * bytes) of this
+                                         * object. Since sometimes
+                                         * the size of objects can
+                                         * not be determined exactly
+                                         * (for example: what is the
+                                         * memory consumption of an
+                                         * STL @p{std::map} type with a
+                                         * certain number of
+                                         * elements?), this is only
+                                         * an estimate. however often
+                                         * quite close to the true
+                                         * value.
+                                         */
+       unsigned int memory_consumption () const;
+    };
+
+                                    /**
+                                     * Flags describing the details
+                                     * of output in VTK format. At
+                                     * present no flags are
+                                     * implemented.
+                                     */
+    struct VtkFlags 
     {
       private:
                                         /**
@@ -1058,6 +1123,19 @@ class DataOutBase
     static void write_gmv (const typename std::vector<Patch<dim,spacedim> > &patches,
                           const std::vector<std::string>          &data_names,
                           const GmvFlags                          &flags,
+                          std::ostream                            &out);
+
+                                    /**
+                                     * Write the given list of
+                                     * patches to the output stream
+                                     * in vtk format. See the general
+                                     * documentation for more
+                                     * information on the parameters.
+                                     */
+    template <int dim, int spacedim>
+    static void write_vtk (const typename std::vector<Patch<dim,spacedim> > &patches,
+                          const std::vector<std::string>          &data_names,
+                          const VtkFlags                          &flags,
                           std::ostream                            &out);
 
                                     /**
@@ -1127,7 +1205,8 @@ class DataOutBase
                                      * coordinates, one height value)
                                      * to the direction of sight.
                                      */
-    class EpsCell2d {
+    class EpsCell2d
+    {
       public:
        
                                         /**
@@ -1172,6 +1251,11 @@ class DataOutBase
                                      * from the main thread, and is
                                      * thus moved into this separate
                                      * function.
+                                     *
+                                     * Note that because of the close
+                                     * similarity of the two formats,
+                                     * this function is also used by
+                                     * the Vtk output function.
                                      */
     template <int dim, int spacedim>
     static void
@@ -1299,96 +1383,122 @@ class DataOutInterface : private DataOutBase
 {
   public:
                                     /**
-                                     * Provide a data type specifying the
-                                     * presently supported output formats.
+                                     * Provide a data type specifying
+                                     * the presently supported output
+                                     * formats.
                                      */
-    enum OutputFormat { default_format, ucd, gnuplot, povray, eps, gmv };
+    enum OutputFormat { default_format, ucd, gnuplot, povray, eps, gmv, vtk };
 
                                     /**
-                                     * Obtain data through the @p{get_patches}
-                                     * function and write it to @p{out} in
-                                     * UCD format.
+                                     * Obtain data through the
+                                     * @p{get_patches} function and
+                                     * write it to @p{out} in UCD
+                                     * format.
                                      */
     void write_ucd (std::ostream &out) const;
 
                                     /**
-                                     * Obtain data through the @p{get_patches}
-                                     * function and write it to @p{out} in
-                                     * GNUPLOT format.
+                                     * Obtain data through the
+                                     * @p{get_patches} function and
+                                     * write it to @p{out} in GNUPLOT
+                                     * format.
                                      */
     void write_gnuplot (std::ostream &out) const;
 
                                     /**
-                                     * Obtain data through the @p{get_patches}
-                                     * function and write it to @p{out} in
-                                     * POVRAY format.
+                                     * Obtain data through the
+                                     * @p{get_patches} function and
+                                     * write it to @p{out} in POVRAY
+                                     * format.
                                      */
     void write_povray (std::ostream &out) const;
 
                                     /**
-                                     * Obtain data through the @p{get_patches}
-                                     * function and write it to @p{out} in
-                                     * EPS format.
+                                     * Obtain data through the
+                                     * @p{get_patches} function and
+                                     * write it to @p{out} in EPS
+                                     * format.
                                      */
     void write_eps (std::ostream &out) const;
 
                                     /**
-                                     * Obtain data through the @p{get_patches}
-                                     * function and write it to @p{out} in
-                                     * GMV format.
+                                     * Obtain data through the
+                                     * @p{get_patches} function and
+                                     * write it to @p{out} in GMV
+                                     * format.
                                      */
     void write_gmv (std::ostream &out) const;
 
+                                    /**
+                                     * Obtain data through the
+                                     * @p{get_patches} function and
+                                     * write it to @p{out} in Vtk
+                                     * format.
+                                     */
+    void write_vtk (std::ostream &out) const;
+
                                     /**
-                                     * Write data and grid to @p{out} according
-                                     * to the given data format. This function
-                                     * simply calls the appropriate
-                                     * @p{write_*} function. If no output format is
-                                     * requested, the @p{default_format} is written.
+                                     * Write data and grid to @p{out}
+                                     * according to the given data
+                                     * format. This function simply
+                                     * calls the appropriate
+                                     * @p{write_*} function. If no
+                                     * output format is requested,
+                                     * the @p{default_format} is
+                                     * written.
                                      *
-                                     * An error occurs if no format is provided and
-                                     * the default format is @p{default_format}.
+                                     * An error occurs if no format
+                                     * is provided and the default
+                                     * format is @p{default_format}.
                                      */
     void write (std::ostream       &out,
                const OutputFormat  output_format = default_format) const;
 
                                     /**
-                                     * Set the default format. The value set here
-                                     * is used anytime, output for format @p{default_format} is
+                                     * Set the default format. The
+                                     * value set here is used
+                                     * anytime, output for format
+                                     * @p{default_format} is
                                      * requested.
                                      */
-    void set_default_format(const OutputFormat default_format);
+    void set_default_format (const OutputFormat default_format);
 
                                     /**
-                                     * Set the flags to be used for output
-                                     * in UCD 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.
+                                     * 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 POVRAY format.
+                                     * Set the flags to be used for
+                                     * output in POVRAY format.
                                      */
     void set_flags (const PovrayFlags &povray_flags);
 
                                     /**
-                                     * Set the flags to be used for output
-                                     * in EPS output.
+                                     * Set the flags to be used for
+                                     * output in EPS output.
                                      */
     void set_flags (const EpsFlags &eps_flags);
 
                                     /**
-                                     * Set the flags to be used for output
-                                     * in GMV format.
+                                     * Set the flags to be used for
+                                     * output in GMV format.
                                      */
     void set_flags (const GmvFlags &gmv_flags);
 
+                                    /**
+                                     * Set the flags to be used for
+                                     * output in GMV format.
+                                     */
+    void set_flags (const VtkFlags &vtk_flags);
+    
 
                                     /**
                                      * Provide a function which tells us which
@@ -1400,82 +1510,97 @@ class DataOutInterface : private DataOutBase
                                      * @item @p{gnuplot}: @p{.gnuplot}
                                      * @item @p{povray}: @p{.pov}
                                      * @item @p{eps}: @p{.eps}
-                                     * @item @p{gmv}: @p{.gmv}.
+                                     * @item @p{gmv}: @p{.gmv}
+                                     * @item @p{vtk}: @p{.vtk}.
                                      * @end{itemize}
                                      *
                                      * If this function is called
-                                     * with no argument or @p{default_format}, the
-                                     * suffix for the
-                                     * @p{default_format} is returned.
+                                     * with no argument or
+                                     * @p{default_format}, the suffix
+                                     * for the @p{default_format} is
+                                     * returned.
                                      */
     std::string default_suffix (const OutputFormat output_format = default_format) const;
 
                                     /**
-                                     * Return the @p{OutputFormat} value
-                                     * corresponding to the given string. If
-                                     * the string does not match any known
-                                     * format, an exception is thrown.
+                                     * Return the @p{OutputFormat}
+                                     * value corresponding to the
+                                     * given string. If the string
+                                     * does not match any known
+                                     * format, an exception is
+                                     * thrown.
                                      *
-                                     * Since this function does not need data
-                                     * from this object, it is static and can
-                                     * thus be called without creating an
-                                     * object of this class. Its main purpose
-                                     * is to allow a program to use any
-                                     * implemented output format without the
-                                     * need to extend the program's parser
-                                     * each time a new format is implemented.
+                                     * Since this function does not
+                                     * need data from this object, it
+                                     * is static and can thus be
+                                     * called without creating an
+                                     * object of this class. Its main
+                                     * purpose is to allow a program
+                                     * to use any implemented output
+                                     * format without the need to
+                                     * extend the program's parser
+                                     * each time a new format is
+                                     * implemented.
                                      *
-                                     * To get a list of presently available
-                                     * format names, e.g. to give it to the
-                                     * @p{ParameterHandler} class, use the
-                                     * function @p{get_output_format_names ()}.
+                                     * To get a list of presently
+                                     * available format names,
+                                     * e.g. to give it to the
+                                     * @p{ParameterHandler} class,
+                                     * use the function
+                                     * @p{get_output_format_names
+                                     * ()}.
                                      */
     static OutputFormat parse_output_format (const std::string &format_name);
 
                                     /**
-                                     * Return a list of implemented output
-                                     * formats. The different names are
-                                     * separated by vertical bar signs (@p{`|'})
-                                     * as used by the @p{ParameterHandler}
-                                     * classes.
+                                     * Return a list of implemented
+                                     * output formats. The different
+                                     * names are separated by
+                                     * vertical bar signs (@p{`|'})
+                                     * as used by the
+                                     * @p{ParameterHandler} classes.
                                      */
     static std::string get_output_format_names ();
 
                                     /**
-                                     * Declare parameters for all output
-                                     * formats by declaring subsections
-                                     * within the parameter file for each
-                                     * output format and call the
-                                     * respective @p{declare_parameters}
-                                     * functions of the flag classes for
-                                     * each output format.
+                                     * Declare parameters for all
+                                     * output formats by declaring
+                                     * subsections within the
+                                     * parameter file for each output
+                                     * format and call the respective
+                                     * @p{declare_parameters}
+                                     * functions of the flag classes
+                                     * for each output format.
                                      *
-                                     * Some of the declared subsections
-                                     * may not contain entries, if the
-                                     * respective format does not
-                                     * export any flags.
+                                     * Some of the declared
+                                     * subsections may not contain
+                                     * entries, if the respective
+                                     * format does not export any
+                                     * flags.
                                      *
-                                     * Note that the top-level parameters
-                                     * denoting the number of subdivisions
-                                     * per patch and the output format
-                                     * are not declared, since they are
-                                     * only passed to virtual functions
-                                     * and are not stored inside objects
-                                     * of this type. You have to
-                                     * declare them yourself.
+                                     * Note that the top-level
+                                     * parameters denoting the number
+                                     * of subdivisions per patch and
+                                     * the output format are not
+                                     * declared, since they are only
+                                     * passed to virtual functions
+                                     * and are not stored inside
+                                     * objects of this type. You have
+                                     * to declare them yourself.
                                      */
     static void declare_parameters (ParameterHandler &prm);
 
                                     /**
-                                     * Read the parameters declared in
-                                     * @p{declare_parameters} and set the
-                                     * flags for the output formats
-                                     * accordingly.
+                                     * Read the parameters declared
+                                     * in @p{declare_parameters} and
+                                     * set the flags for the output
+                                     * formats accordingly.
                                      *
-                                     * The flags thus obtained overwrite
-                                     * all previous contents of the flag
-                                     * objects as default-constructed or
-                                     * set by the @p{set_flags} function.
+                                     * The flags thus obtained
+                                     * overwrite all previous
+                                     * contents of the flag objects
+                                     * as default-constructed or set
+                                     * by the @p{set_flags} function.
                                      */
     void parse_parameters (ParameterHandler &prm);
     
@@ -1504,24 +1629,27 @@ class DataOutInterface : private DataOutBase
     
   protected:
                                     /**
-                                     * This is the abstract function through
-                                     * which derived classes propagate
-                                     * preprocessed data in the form of
-                                     * @p{Patch} structures (declared in
-                                     * the base class @p{DataOutBase}) to
-                                     * the actual output function. You
-                                     * need to overload this function to
-                                     * allow the output functions to
-                                     * know what they shall print.
+                                     * This is the abstract function
+                                     * through which derived classes
+                                     * propagate preprocessed data in
+                                     * the form of @p{Patch}
+                                     * structures (declared in the
+                                     * base class @p{DataOutBase}) to
+                                     * the actual output
+                                     * function. You need to overload
+                                     * this function to allow the
+                                     * output functions to know what
+                                     * they shall print.
                                      */
     virtual const typename std::vector<typename DataOutBase::Patch<dim,spacedim> > &
     get_patches () const = 0;
 
                                     /**
-                                     * Abstract virtual function through
-                                     * which the names of data sets are
-                                     * obtained by the output functions
-                                     * of the base class.
+                                     * Abstract virtual function
+                                     * through which the names of
+                                     * data sets are obtained by the
+                                     * output functions of the base
+                                     * class.
                                      */
     virtual std::vector<std::string> get_dataset_names () const = 0;
 
@@ -1537,41 +1665,55 @@ class DataOutInterface : private DataOutBase
     OutputFormat default_fmt;
     
                                     /**
-                                     * Flags to be used upon output of UCD
-                                     * data. Can be changed by using the
-                                     * @p{set_flags} function.
+                                     * Flags to be used upon output
+                                     * of UCD data. Can be changed by
+                                     * using the @p{set_flags}
+                                     * function.
                                      */
     UcdFlags     ucd_flags;
 
                                     /**
-                                     * Flags to be used upon output of GNUPLOT
-                                     * data. Can be changed by using the
+                                     * Flags to be used upon output
+                                     * of GNUPLOT data. Can be
+                                     * changed by using the
                                      * @p{set_flags} function.
                                      */
     GnuplotFlags gnuplot_flags;
 
                                     /**
-                                     * Flags to be used upon output of POVRAY
-                                     * data. Can be changed by using the
-                                     * @p{set_flags} function.
+                                     * Flags to be used upon output
+                                     * of POVRAY data. Can be changed
+                                     * by using the @p{set_flags}
+                                     * function.
                                      */
     PovrayFlags povray_flags;
 
                                     /**
-                                     * Flags to be used upon output of EPS
-                                     * data in one space dimension. Can be
-                                     * changed by using the @p{set_flags}
+                                     * Flags to be used upon output
+                                     * of EPS data in one space
+                                     * dimension. Can be changed by
+                                     * using the @p{set_flags}
                                      * function.
                                      */
     EpsFlags     eps_flags;    
 
                                     /**
-                                     * Flags to be used upon output of gmv
-                                     * data in one space dimension. Can be
-                                     * changed by using the @p{set_flags}
+                                     * Flags to be used upon output
+                                     * of gmv data in one space
+                                     * dimension. Can be changed by
+                                     * using the @p{set_flags}
                                      * function.
                                      */
     GmvFlags     gmv_flags;
+
+                                    /**
+                                     * Flags to be used upon output
+                                     * of vtk data in one space
+                                     * dimension. Can be changed by
+                                     * using the @p{set_flags}
+                                     * function.
+                                     */
+    VtkFlags     vtk_flags;
 };
 
 
index 618cb56fc12a8b237ea39344aa9e822e4aef613f..395f2970feae754c6285444215f30385f01626a3 100644 (file)
@@ -327,6 +327,25 @@ DataOutBase::GmvFlags::memory_consumption () const
 
 
 
+void DataOutBase::VtkFlags::declare_parameters (ParameterHandler &/*prm*/)
+{};
+
+
+
+void DataOutBase::VtkFlags::parse_parameters (ParameterHandler &/*prm*/)
+{};
+
+
+unsigned int
+DataOutBase::VtkFlags::memory_consumption () const
+{
+                                  // only simple data elements, so
+                                  // use sizeof operator
+  return sizeof (*this);
+};
+
+
+
 unsigned int DataOutBase::memory_consumption ()
 {
   return 0;
index 7e85c1bc5ddb27ee6053c55dd502e403314889d5..796b2c6d1f4151c62ce33924cc6949a0e9b0140e 100644 (file)
@@ -1768,6 +1768,355 @@ void DataOutBase::write_gmv (const typename std::vector<Patch<dim,spacedim> > &p
 
 
 
+template <int dim, int spacedim>
+void DataOutBase::write_vtk (const typename std::vector<Patch<dim,spacedim> > &patches,
+                            const std::vector<std::string>          &data_names,
+                            const VtkFlags                          &/*flags*/,
+                            std::ostream                            &out) 
+{
+  AssertThrow (out, ExcIO());
+
+  Assert (patches.size() > 0, ExcNoPatches());
+  const unsigned int n_data_sets = data_names.size();
+                                  // check against # of data sets in
+                                  // first patch. checks against all
+                                  // other patches are made in
+                                  // write_gmv_reorder_data_vectors
+  Assert (n_data_sets == patches[0].data.m(),
+         ExcUnexpectedNumberOfDatasets (patches[0].data.m(), n_data_sets));
+  
+  
+                                  ///////////////////////
+                                  // preamble
+  if (true)
+    {
+      time_t  time1= time (0);
+      tm     *time = localtime(&time1);
+      out << "# vtk DataFile Version 3.0"
+         << std::endl
+         << "This file was generated by the deal.II library on "
+         << time->tm_year+1900 << "/"
+         << time->tm_mon+1 << "/"
+         << time->tm_mday << " at "
+         << time->tm_hour << ":"
+         << std::setw(2) << time->tm_min << ":"
+         << std::setw(2) << time->tm_sec
+         << std::endl
+         << "ASCII"
+         << std::endl
+         << "DATASET UNSTRUCTURED_GRID"
+         << std::endl;
+    };
+  
+
+                                  // first count the number of cells
+                                  // and cells for later use
+  unsigned int n_cells = 0,
+              n_nodes = 0;
+  for (typename std::vector<Patch<dim,spacedim> >::const_iterator patch=patches.begin();
+       patch!=patches.end(); ++patch)
+    switch (dim)
+      {
+       case 1:
+             n_cells += patch->n_subdivisions;
+             n_nodes += patch->n_subdivisions+1;
+             break;
+       case 2:
+             n_cells += patch->n_subdivisions *
+                        patch->n_subdivisions;
+             n_nodes += (patch->n_subdivisions+1) *
+                        (patch->n_subdivisions+1);
+             break;
+       case 3:
+             n_cells += patch->n_subdivisions *
+                        patch->n_subdivisions *
+                        patch->n_subdivisions;
+             n_nodes += (patch->n_subdivisions+1) *
+                        (patch->n_subdivisions+1) *
+                        (patch->n_subdivisions+1);
+             break;
+       default:
+             Assert (false, ExcNotImplemented());
+      };
+
+
+                                  // in gmv format the vertex
+                                  // coordinates and the data have an
+                                  // order that is a bit unpleasant
+                                  // (first all x coordinates, then
+                                  // all y coordinate, ...; first all
+                                  // data of variable 1, then
+                                  // variable 2, etc), so we have to
+                                  // copy the data vectors a bit around
+                                  //
+                                  // note that we copy vectors when
+                                  // looping over the patches since we
+                                  // have to write them one variable
+                                  // at a time and don't want to use
+                                  // more than one loop
+                                  //
+                                  // this copying of data vectors can
+                                  // be done while we already output
+                                  // the vertices, so do this on a
+                                  // separate thread and when wanting
+                                  // to write out the data, we wait
+                                  // for that thread to finish
+  std::vector<std::vector<double> > data_vectors (n_data_sets,
+                                                 std::vector<double> (n_nodes));
+  Threads::ThreadManager thread_manager;
+  Threads::spawn (thread_manager,
+                 Threads::encapsulate (&DataOutBase::template
+                                       write_gmv_reorder_data_vectors<dim,spacedim>)
+                 .collect_args(patches, data_vectors));
+
+                                  ///////////////////////////////
+                                  // first make up a list of used
+                                  // vertices along with their
+                                  // coordinates
+                                  //
+                                  // note that we have to print
+                                  // d=1..3 dimensions
+  out << "POINTS " << n_nodes << " double" << std::endl;
+  for (unsigned int d=1; d<=3; ++d)
+    {
+      for (typename std::vector<Patch<dim,spacedim> >::const_iterator patch=patches.begin();
+          patch!=patches.end(); ++patch)
+       {
+         const unsigned int n_subdivisions = patch->n_subdivisions;
+         
+                                          // if we have nonzero values for
+                                          // this coordinate
+         if (d<=spacedim)
+           {
+             switch (dim)
+               {
+                 case 1:
+                 {
+                   for (unsigned int i=0; i<n_subdivisions+1; ++i)
+                     out << ((patch->vertices[1](0) * i / n_subdivisions) +
+                             (patch->vertices[0](0) * (n_subdivisions-i) / n_subdivisions))
+                         << ' ';
+                   break;
+                 };
+                  
+                 case 2:
+                 {
+                   for (unsigned int i=0; i<n_subdivisions+1; ++i)
+                     for (unsigned int j=0; j<n_subdivisions+1; ++j)
+                       {
+                         const double x_frac = i * 1./n_subdivisions,
+                                      y_frac = j * 1./n_subdivisions;
+                         
+                                                          // compute coordinates for
+                                                          // this patch point
+                         out << (((patch->vertices[1](d-1) * x_frac) +
+                                  (patch->vertices[0](d-1) * (1-x_frac))) * (1-y_frac) +
+                                 ((patch->vertices[2](d-1) * x_frac) +
+                                  (patch->vertices[3](d-1) * (1-x_frac))) * y_frac)
+                             << ' ';
+                       };
+                   break;
+                 };
+                  
+                 case 3:
+                 {
+                   for (unsigned int i=0; i<n_subdivisions+1; ++i)
+                     for (unsigned int j=0; j<n_subdivisions+1; ++j)
+                       for (unsigned int k=0; k<n_subdivisions+1; ++k)
+                         {
+                                                            // note the broken
+                                                            // design of hexahedra
+                                                            // in deal.II, where
+                                                            // first the z-component
+                                                            // is counted up, before
+                                                            // increasing the y-
+                                                            // coordinate
+                           const double x_frac = i * 1./n_subdivisions,
+                                        y_frac = k * 1./n_subdivisions,
+                                        z_frac = j * 1./n_subdivisions;
+                           
+                                                            // compute coordinates for
+                                                            // this patch point
+                           out << ((((patch->vertices[1](d-1) * x_frac) +
+                                     (patch->vertices[0](d-1) * (1-x_frac))) * (1-y_frac) +
+                                    ((patch->vertices[2](d-1) * x_frac) +
+                                     (patch->vertices[3](d-1) * (1-x_frac))) * y_frac)   * (1-z_frac) +
+                                   (((patch->vertices[5](d-1) * x_frac) +
+                                     (patch->vertices[4](d-1) * (1-x_frac))) * (1-y_frac) +
+                                    ((patch->vertices[6](d-1) * x_frac) +
+                                     (patch->vertices[7](d-1) * (1-x_frac))) * y_frac)   * z_frac)
+                               << ' ';
+                         };
+             
+                   break;
+                 };
+                  
+                 default:
+                       Assert (false, ExcNotImplemented());
+               };
+           }
+         else
+                                            // d>spacedim. write zeros instead
+           {
+             const unsigned int n_points
+               = static_cast<unsigned int>(pow (static_cast<double>(n_subdivisions+1), dim));
+             for (unsigned int i=0; i<n_points; ++i)
+               out << "0 ";
+           };
+       };
+      out << std::endl;
+    };
+
+  out << std::endl;
+
+                                  /////////////////////////////////
+                                  // now for the cells. note that
+                                  // vertices are counted from 1 onwards
+  if (true)
+    {
+      out << "CELLS " << n_cells << ' '
+         << n_cells*(GeometryInfo<dim>::vertices_per_cell+1)
+         << std::endl;
+
+
+      unsigned int first_vertex_of_patch = 0;
+      
+      for (typename std::vector<Patch<dim,spacedim> >::const_iterator patch=patches.begin();
+          patch!=patches.end(); ++patch)
+       {
+         const unsigned int n_subdivisions = patch->n_subdivisions;
+
+                                          // write out the cells making
+                                          // up this patch
+         switch (dim)
+           {
+             case 1:
+             {
+               for (unsigned int i=0; i<n_subdivisions; ++i)
+                 out << "2 "
+                     << first_vertex_of_patch+i+1 << ' '
+                     << first_vertex_of_patch+i+1+1 << std::endl;
+               break;
+             };
+              
+             case 2:
+             {
+               for (unsigned int i=0; i<n_subdivisions; ++i)
+                 for (unsigned int j=0; j<n_subdivisions; ++j)
+                   out << "4 "
+                       << first_vertex_of_patch+i*(n_subdivisions+1)+j+1 << ' '
+                       << first_vertex_of_patch+(i+1)*(n_subdivisions+1)+j+1 << ' '
+                       << first_vertex_of_patch+(i+1)*(n_subdivisions+1)+j+1+1 << ' '
+                       << first_vertex_of_patch+i*(n_subdivisions+1)+j+1+1
+                       << std::endl;
+               break;
+             };
+              
+             case 3:
+             {
+               for (unsigned int i=0; i<n_subdivisions; ++i)
+                 for (unsigned int j=0; j<n_subdivisions; ++j)
+                   for (unsigned int k=0; k<n_subdivisions; ++k)
+                     {
+                       out << "8 "
+                                                          // note: vertex indices start with 1!
+                           << first_vertex_of_patch+(i*(n_subdivisions+1)+j      )*(n_subdivisions+1)+k  +1 << ' '
+                           << first_vertex_of_patch+((i+1)*(n_subdivisions+1)+j  )*(n_subdivisions+1)+k  +1 << ' '
+                           << first_vertex_of_patch+((i+1)*(n_subdivisions+1)+j+1)*(n_subdivisions+1)+k  +1 << ' '
+                           << first_vertex_of_patch+(i*(n_subdivisions+1)+j+1    )*(n_subdivisions+1)+k  +1 << ' '
+                           << first_vertex_of_patch+(i*(n_subdivisions+1)+j      )*(n_subdivisions+1)+k+1+1 << ' '
+                           << first_vertex_of_patch+((i+1)*(n_subdivisions+1)+j  )*(n_subdivisions+1)+k+1+1 << ' '
+                           << first_vertex_of_patch+((i+1)*(n_subdivisions+1)+j+1)*(n_subdivisions+1)+k+1+1 << ' '
+                           << first_vertex_of_patch+(i*(n_subdivisions+1)+j+1    )*(n_subdivisions+1)+k+1+1 << ' '
+                           << std::endl;
+                     };
+               break;
+             };
+
+             default:
+                   Assert (false, ExcNotImplemented());
+           };
+
+
+                                          // finally update the number
+                                          // of the first vertex of this patch
+         switch (dim)
+           {
+             case 1:
+                   first_vertex_of_patch += n_subdivisions+1;
+                   break;
+             case 2:
+                   first_vertex_of_patch += (n_subdivisions+1) *
+                                            (n_subdivisions+1);
+                   break;
+             case 3:
+                   first_vertex_of_patch += (n_subdivisions+1) *
+                                            (n_subdivisions+1) *
+                                            (n_subdivisions+1);
+                   break;
+             default:
+                   Assert (false, ExcNotImplemented());
+           };
+       };
+      out << std::endl;
+
+                                      // next output the types of the
+                                      // cells. since all cells are
+                                      // the same, this is simple
+      for (unsigned int i=0; i<n_cells; ++i)
+       switch (dim)
+         {
+           case 1:
+                 out << "3\n";    // represents VTK_LINE
+                 break;
+           case 2:
+                 out << "9\n";    // represents VTK_QUAD
+                 break;
+           case 3:
+                 out << "12\n";    // represents VTK_HEXAHEDRON
+                 break;
+           default:
+                 Assert (false, ExcNotImplemented());
+         };
+    };
+
+                                  ///////////////////////////////////////
+                                  // data output.
+
+                                  // now write the data vectors to
+                                  // @p{out} first make sure that all
+                                  // data is in place
+  thread_manager.wait ();
+
+                                  // then write data.
+                                  // the '1' means: node data (as opposed
+                                  // to cell data, which we do not
+                                  // support explicitely here)
+  for (unsigned int data_set=0; data_set<n_data_sets; ++data_set)
+    {
+      out << "POINTDATA " << n_nodes
+         << std::endl
+         << "SCALARS "
+         << data_names[data_set]
+         << " double 1"
+         << std::endl
+         << "LOOKUP_TABLE default"
+         << std::endl;
+      copy(data_vectors[data_set].begin(),
+          data_vectors[data_set].end(),
+          std::ostream_iterator<double>(out, " "));
+      out << std::endl
+         << std::endl;
+    };
+  
+                                  // assert the stream is still ok
+  AssertThrow (out, ExcIO());
+};
+
+
+
+
 template <int dim, int spacedim>
 void
 DataOutBase::write_gmv_reorder_data_vectors (const typename std::vector<Patch<dim,spacedim> > &patches,
@@ -1904,6 +2253,15 @@ void DataOutInterface<dim,spacedim>::write_gmv (std::ostream &out) const
 
 
 
+template <int dim, int spacedim>
+void DataOutInterface<dim,spacedim>::write_vtk (std::ostream &out) const 
+{
+  DataOutBase::write_vtk (get_patches(), get_dataset_names(),
+                         vtk_flags, out);
+};
+
+
+
 template <int dim, int spacedim>
 void DataOutInterface<dim,spacedim>::write (std::ostream &out,
                                            OutputFormat output_format) const
@@ -1933,6 +2291,10 @@ void DataOutInterface<dim,spacedim>::write (std::ostream &out,
            write_gmv (out);
            break;
            
+      case vtk:
+           write_vtk (out);
+           break;
+           
       default:
            Assert (false, ExcNotImplemented());
     };
@@ -1989,6 +2351,14 @@ void DataOutInterface<dim,spacedim>::set_flags (const GmvFlags &flags)
 
 
 
+template <int dim, int spacedim>
+void DataOutInterface<dim,spacedim>::set_flags (const VtkFlags &flags) 
+{
+  vtk_flags = flags;
+};
+
+
+
 template <int dim, int spacedim>
 std::string DataOutInterface<dim,spacedim>::default_suffix (OutputFormat output_format) const
 {
@@ -2012,6 +2382,9 @@ std::string DataOutInterface<dim,spacedim>::default_suffix (OutputFormat output_
       case gmv:
            return ".gmv";
            
+      case vtk:
+           return ".vtk";
+           
       default: 
            Assert (false, ExcNotImplemented()); 
            return "";
@@ -2039,6 +2412,9 @@ DataOutInterface<dim,spacedim>::parse_output_format (const std::string &format_n
   if (format_name == "gmv")
     return gmv;
   
+  if (format_name == "vtk")
+    return vtk;
+  
   AssertThrow (false, ExcInvalidState ());
 
                                   // return something invalid
@@ -2050,7 +2426,7 @@ DataOutInterface<dim,spacedim>::parse_output_format (const std::string &format_n
 template <int dim, int spacedim>
 std::string DataOutInterface<dim,spacedim>::get_output_format_names ()
 {
-  return "ucd|gnuplot|povray|eps|gmv";
+  return "ucd|gnuplot|povray|eps|gmv|vtk";
 }
 
 
@@ -2080,6 +2456,10 @@ void DataOutInterface<dim,spacedim>::declare_parameters (ParameterHandler &prm)
   prm.enter_subsection ("Gmv output parameters");
   GmvFlags::declare_parameters (prm);
   prm.leave_subsection ();
+
+  prm.enter_subsection ("Vtk output parameters");
+  VtkFlags::declare_parameters (prm);
+  prm.leave_subsection ();
 }
 
 
@@ -2109,6 +2489,10 @@ void DataOutInterface<dim,spacedim>::parse_parameters (ParameterHandler &prm)
   prm.enter_subsection ("Gmv output parameters");
   gmv_flags.parse_parameters (prm);
   prm.leave_subsection ();
+
+  prm.enter_subsection ("Vtk output parameters");
+  vtk_flags.parse_parameters (prm);
+  prm.leave_subsection ();
 }
 
 
@@ -2121,7 +2505,8 @@ DataOutInterface<dim,spacedim>::memory_consumption () const
          MemoryConsumption::memory_consumption (gnuplot_flags) +
          MemoryConsumption::memory_consumption (povray_flags) +
          MemoryConsumption::memory_consumption (eps_flags) +
-         MemoryConsumption::memory_consumption (gmv_flags));
+         MemoryConsumption::memory_consumption (gmv_flags) +
+         MemoryConsumption::memory_consumption (vtk_flags));
 };
 
 
index 22b71f0effbd6505f225a509f787621bfb611ebc..190cc19d063af9730c77c9bf45efae652536ed7a 100644 (file)
@@ -119,6 +119,15 @@ documentation, etc</a>.
 <h3>base</h3>
 
 <ol>
+  <li> <p>
+       We now support the output format of the <a
+       href="http://www.kitware.com/vtk.html"
+       target="_top">Visualization Toolkit (Vtk)</a> from the <code
+       class="class">DataOutBase</code> class and all derived classes.
+       <br>
+       (WB 2001/03/19)
+       </p>
+
   <li> <p>
        New: The class <code
        class="class">TensorProductPolynomials&lt;dim&gt;</code>
index cf64eb8cd87ce8a473278fc1a4c6ddd48696771c..170e67e401e1792ec558714a34268e64ec8aab4e 100644 (file)
            GNUPLOT, 
           <a href="http://www-xdiv.lanl.gov/XCM/gmv/" target="_top">GMV
                                      (general mesh viewer)</a>, 
+           <a href="http://www.kitware.com/vtk.html"
+           target="_top">Visualization Toolkit (Vtk)</a>,
           <a href="http://www.avs.com" target="_top">AVS Explorer</a>,
           <a href="http://www.povray.org" target="_top">Povray</a>,
            and directly to Encapsulated Postscript. GNUPLOT and a
            Unix flavors.
           </p>
     </ul>
-    <p>
 
+    <p>
     In case you didn't find your favorite graphics format above it
     might be of interest that it is quite simple to add it to the list
     of supported ones (take a look at
     base/include/base/data_out_base.h and
     base/source/data_out_base.cc to see how output formats are
-    implemented).
+    implemented), as only a simple intermediate format needs to be
+    converted into a graphics format, without references to cells,
+    nodes, types of finite elements, and the like.
+    </p>
 
 
 

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.