]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add deal.II intermediate form.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 5 Apr 2005 19:18:36 +0000 (19:18 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Tue, 5 Apr 2005 19:18:36 +0000 (19:18 +0000)
git-svn-id: https://svn.dealii.org/trunk@10357 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/include/base/data_out_base.h
deal.II/base/source/data_out_base.cc

index e556fcf7e3844b37befc1444c806be8a3ae90700..127702b7c991346fda64f92b0960c93049df407e 100644 (file)
@@ -324,6 +324,39 @@ class ParameterHandler;
  * there for more information.
  *
  *
+ * @subsection DataOutBaseD2 deal.II intermediate format
+ *
+ * In addition to all the other formats, this class can also write
+ * data in the deal.II format. This is not a format that is read by
+ * any other graphics program, but is rather a direct dump of the
+ * intermediate internal format used by deal.II. This internal format
+ * is generated by the various classes that can generate output using
+ * the DataOutBase class, for example from a finite element solution,
+ * and is then converted in the present class to the final graphics
+ * format. The reason why we offer to write out this intermediate
+ * format is that it can be read back into a deal.II program, which is
+ * helpful in at least two contexts: First, this can be used to later
+ * generate graphical output in any other graphics format presently
+ * understood; this way, it is not necessary to know at run-time which
+ * output format is requested, or if multiple output files in
+ * different formats are needed. Secondly, in contrast to almost all
+ * other graphics formats, it is possible to merge several files that
+ * contain intermediate format data, and generate a single output file
+ * from it, which may be again in intermediate format or any of the
+ * final formats. This latter option is most helpful for parallel
+ * programs: as demonstrated in the step-17 example program, it is
+ * possible to let only one processor generate the graphical output
+ * for the entire parallel program, but this can become vastly
+ * inefficient if many processors are involved, because the load is no
+ * longer balanced. The way out is to let each processor generate
+ * intermediate graphical output for its chunk of the domain, and the
+ * later merge the different files into one, which is an operation
+ * that is much cheaper than the generation of the intermediate data.
+ *
+ * Intermediate format deal.II data is usually stored in files with
+ * the ending <tt>.d2</tt>.
+ *
+ *
  * @section DataOutBaseOP Output parameters
  *
  * All functions take a parameter which is a structure of type <tt>XFlags</tt>, where
@@ -518,7 +551,9 @@ class DataOutBase
        unsigned int memory_consumption () const;
        
                                         /**
-                                         * Value for no neighbor.
+                                         * Value to be used if this
+                                         * patch has no neighbor on
+                                         * one side.
                                          */
        static const unsigned int no_neighbor = deal_II_numbers::invalid_unsigned_int;
                                         /** @addtogroup Exceptions
@@ -1290,6 +1325,66 @@ class DataOutBase
        unsigned int memory_consumption () const;
     };
 
+
+
+                                    /**
+                                     * Flags describing the details
+                                     * of output in deal.II
+                                     * intermediate format. At
+                                     * present no flags are
+                                     * implemented.
+                                     */
+    struct Deal_II_IntermediateFlags
+    {
+      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
+                                         * <tt>private</tt> 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
+                                         * <tt>declare_parameters</tt> 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 <tt>std::map</tt> 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;
+    };
+    
                                     /**
                                      * Write the given list of patches
                                      * to the output stream in OpenDX
@@ -1415,6 +1510,20 @@ class DataOutBase
                           const VtkFlags                          &flags,
                           std::ostream                            &out);
 
+                                    /**
+                                     * Write the given list of
+                                     * patches to the output stream
+                                     * in deal.II intermediate
+                                     * format. See the general
+                                     * documentation for more
+                                     * information on the parameters.
+                                     */
+    template <int dim, int spacedim>
+    static void write_deal_II_intermediate (const std::vector<Patch<dim,spacedim> > &patches,
+                                           const std::vector<std::string>          &data_names,
+                                           const Deal_II_IntermediateFlags         &flags,
+                                           std::ostream                            &out);
+
                                     /**
                                      * Determine an estimate for
                                      * the memory consumption (in
@@ -1699,7 +1808,13 @@ class DataOutInterface : private DataOutBase
                                           /**
                                            * Output in VTK format.
                                            */
-         vtk };
+         vtk,
+                                          /**
+                                           * Output in deal.II
+                                           * intermediate format.
+                                           */
+         deal_II_intermediate
+    };
 
                                      /**
                                       * Destructor. Does nothing, but is
@@ -1783,6 +1898,15 @@ class DataOutInterface : private DataOutBase
                                      */
     void write_vtk (std::ostream &out) const;
 
+                                    /**
+                                     * Obtain data through the
+                                     * <tt>get_patches</tt> function
+                                     * and write it to <tt>out</tt>
+                                     * in deal.II intermediate
+                                     * format.
+                                     */
+    void write_deal_II_intermediate (std::ostream &out) const;
+    
                                     /**
                                      * Write data and grid to <tt>out</tt>
                                      * according to the given data
@@ -1853,9 +1977,15 @@ class DataOutInterface : private DataOutBase
 
                                     /**
                                      * Set the flags to be used for
-                                     * output in GMV format.
+                                     * output in VTK format.
                                      */
     void set_flags (const VtkFlags &vtk_flags);
+
+                                    /**
+                                     * Set the flags to be used for
+                                     * output in VTK format.
+                                     */
+    void set_flags (const Deal_II_IntermediateFlags &deal_II_intermediate_flags);
     
 
                                     /**
@@ -1871,6 +2001,7 @@ class DataOutInterface : private DataOutBase
                                      * <li> <tt>eps</tt>: <tt>.eps</tt>
                                      * <li> <tt>gmv</tt>: <tt>.gmv</tt>
                                      * <li> <tt>vtk</tt>: <tt>.vtk</tt>.
+                                     * <li> <tt>deal_II_intermediate</tt>: <tt>.d2</tt>.
                                      * </ul>
                                      *
                                      * If this function is called
@@ -2084,6 +2215,15 @@ class DataOutInterface : private DataOutBase
                                      * function.
                                      */
     VtkFlags     vtk_flags;
+
+                                    /**
+                                     * Flags to be used upon output
+                                     * of vtk data in one space
+                                     * dimension. Can be changed by
+                                     * using the <tt>set_flags</tt>
+                                     * function.
+                                     */
+    Deal_II_IntermediateFlags     deal_II_intermediate_flags;
 };
 
 
@@ -2098,5 +2238,22 @@ DataOutBase::EpsFlags::RgbValues::is_grey () const
 }
 
 
+/* -------------------- template functions ------------------- */
+
+/**
+ * Output operator for an object of type
+ * <tt>DataOutBase::Patch</tt>. This operator dumps the intermediate
+ * graphics format represented by the patch data structure. It may
+ * later be converted into regular formats for a number of graphics
+ * programs.
+ *
+ * @author Wolfgang Bangerth, 2005
+ */
+template <int dim, int spacedim>
+std::ostream &
+operator << (std::ostream                           &out,
+            const DataOutBase::Patch<dim,spacedim> &patch);
+
+
 
 #endif
index 158d65a77ffa7e3b57ae7abbcac5221ab6a888da..0f940c4ad8179148b58b8efac2a0035167d9312b 100644 (file)
@@ -476,6 +476,25 @@ DataOutBase::VtkFlags::memory_consumption () const
 
 
 
+void DataOutBase::Deal_II_IntermediateFlags::declare_parameters (ParameterHandler &/*prm*/)
+{}
+
+
+
+void DataOutBase::Deal_II_IntermediateFlags::parse_parameters (ParameterHandler &/*prm*/)
+{}
+
+
+unsigned int
+DataOutBase::Deal_II_IntermediateFlags::memory_consumption () const
+{
+                                  // only simple data elements, so
+                                  // use sizeof operator
+  return sizeof (*this);
+}
+
+
+
 unsigned int DataOutBase::memory_consumption ()
 {
   return 0;
@@ -3874,6 +3893,30 @@ void DataOutBase::write_vtk (const std::vector<Patch<dim,spacedim> > &patches,
 
 
 
+template <int dim, int spacedim>
+void
+DataOutBase::
+write_deal_II_intermediate (const std::vector<Patch<dim,spacedim> > &patches,
+                           const std::vector<std::string>          &data_names,
+                           const Deal_II_IntermediateFlags         &/*flags*/,
+                           std::ostream                            &out) 
+{
+  out << "[deal.II intermediate format graphics data]" << std::endl
+      << "[written by deal.II version "
+      << DEAL_II_MAJOR << '.' << DEAL_II_MINOR << "]" << std::endl;
+
+  out << data_names.size() << std::endl;
+  for (unsigned int i=0; i<data_names.size(); ++i)
+    out << data_names[i] << std::endl;
+  
+  out << patches.size() << std::endl;
+  for (unsigned int i=0; i<patches.size(); ++i)
+    out << patches[i] << std::endl;
+
+  out << std::endl;
+} 
+
+
 
 template <int dim, int spacedim>
 void
@@ -4054,6 +4097,16 @@ void DataOutInterface<dim,spacedim>::write_vtk (std::ostream &out) const
 
 
 
+template <int dim, int spacedim>
+void DataOutInterface<dim,spacedim>::
+write_deal_II_intermediate (std::ostream &out) const 
+{
+  DataOutBase::write_deal_II_intermediate (get_patches(), get_dataset_names(),
+                                          deal_II_intermediate_flags, out);
+}
+
+
+
 template <int dim, int spacedim>
 void
 DataOutInterface<dim,spacedim>::write (std::ostream &out,
@@ -4100,6 +4153,10 @@ DataOutInterface<dim,spacedim>::write (std::ostream &out,
       case vtk:
            write_vtk (out);
            break;
+
+      case deal_II_intermediate:
+           write_deal_II_intermediate (out);
+           break;
            
       default:
            Assert (false, ExcNotImplemented());
@@ -4190,6 +4247,15 @@ DataOutInterface<dim,spacedim>::set_flags (const VtkFlags &flags)
 
 
 
+template <int dim, int spacedim>
+void
+DataOutInterface<dim,spacedim>::set_flags (const Deal_II_IntermediateFlags &flags) 
+{
+  deal_II_intermediate_flags = flags;
+}
+
+
+
 template <int dim, int spacedim>
 std::string
 DataOutInterface<dim,spacedim>::
@@ -4228,6 +4294,9 @@ default_suffix (const OutputFormat output_format_) const
       case vtk:
            return ".vtk";
            
+      case deal_II_intermediate:
+           return ".d2";
+           
       default: 
            Assert (false, ExcNotImplemented()); 
            return "";
@@ -4268,6 +4337,9 @@ parse_output_format (const std::string &format_name)
   if (format_name == "vtk")
     return vtk;
   
+  if (format_name == "deal.II intermediate")
+    return deal_II_intermediate;
+  
   AssertThrow (false, ExcInvalidState ());
 
                                   // return something invalid
@@ -4280,7 +4352,7 @@ template <int dim, int spacedim>
 std::string
 DataOutInterface<dim,spacedim>::get_output_format_names ()
 {
-  return "dx|ucd|gnuplot|povray|eps|gmv|tecplot|vtk";
+  return "dx|ucd|gnuplot|povray|eps|gmv|tecplot|vtk|deal.II intermediate";
 }
 
 
@@ -4323,6 +4395,11 @@ DataOutInterface<dim,spacedim>::declare_parameters (ParameterHandler &prm)
   prm.enter_subsection ("Vtk output parameters");
   VtkFlags::declare_parameters (prm);
   prm.leave_subsection ();
+
+
+  prm.enter_subsection ("deal.II intermediate output parameters");
+  Deal_II_IntermediateFlags::declare_parameters (prm);
+  prm.leave_subsection ();
 }
 
 
@@ -4365,6 +4442,10 @@ DataOutInterface<dim,spacedim>::parse_parameters (ParameterHandler &prm)
   prm.enter_subsection ("Vtk output parameters");
   vtk_flags.parse_parameters (prm);
   prm.leave_subsection ();
+
+  prm.enter_subsection ("deal.II intermediate output parameters");
+  deal_II_intermediate_flags.parse_parameters (prm);
+  prm.leave_subsection ();
 }
 
 
@@ -4381,7 +4462,42 @@ DataOutInterface<dim,spacedim>::memory_consumption () const
          MemoryConsumption::memory_consumption (eps_flags) +
          MemoryConsumption::memory_consumption (gmv_flags) +
          MemoryConsumption::memory_consumption (tecplot_flags) +
-         MemoryConsumption::memory_consumption (vtk_flags));
+         MemoryConsumption::memory_consumption (vtk_flags) +
+         MemoryConsumption::memory_consumption (deal_II_intermediate_flags));
+}
+
+
+
+template <int dim, int spacedim>
+std::ostream &
+operator << (std::ostream                           &out,
+            const DataOutBase::Patch<dim,spacedim> &patch)
+{
+                                  // write a header line
+  out << "[deal.II intermediate Patch<" << dim << ',' << spacedim << ">]"
+      << std::endl;
+
+                                  // then write all the data that is
+                                  // in this patch
+  for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
+    out << patch.vertices[i] << " ";
+  out << std::endl;
+
+  for (unsigned int i=0; i<GeometryInfo<dim>::faces_per_cell; ++i)
+    out << patch.neighbors[i] << " ";
+  out << std::endl;
+
+  out << patch.patch_index << ' ' << patch.n_subdivisions
+      << std::endl;
+
+  out << patch.data.n_rows() << ' ' << patch.data.n_cols() << std::endl;
+  for (unsigned int i=0; i<patch.data.n_rows(); ++i)
+    for (unsigned int j=0; j<patch.data.n_cols(); ++j)
+      out << patch.data[i][j] << ' ';
+  out << std::endl;
+  out << std::endl;
+
+  return out;
 }
 
 
@@ -4410,3 +4526,34 @@ template class DataOutBase::Patch<2,3>;
 template class DataOutInterface<3,4>;
 template class DataOutBase::Patch<3,4>;
 
+// output operators
+template
+std::ostream &
+operator << (std::ostream                  &out,
+            const DataOutBase::Patch<1,1> &patch);
+template
+std::ostream &
+operator << (std::ostream                  &out,
+            const DataOutBase::Patch<2,2> &patch);
+template
+std::ostream &
+operator << (std::ostream                  &out,
+            const DataOutBase::Patch<3,3> &patch);
+template
+std::ostream &
+operator << (std::ostream                  &out,
+            const DataOutBase::Patch<4,4> &patch);
+template
+std::ostream &
+operator << (std::ostream                  &out,
+            const DataOutBase::Patch<1,2> &patch);
+template
+std::ostream &
+operator << (std::ostream                  &out,
+            const DataOutBase::Patch<2,3> &patch);
+template
+std::ostream &
+operator << (std::ostream                  &out,
+            const DataOutBase::Patch<3,4> &patch);
+
+

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.