]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add the possibility to identify certain output components as vectors. Still has to...
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Sun, 14 Oct 2007 05:25:26 +0000 (05:25 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Sun, 14 Oct 2007 05:25:26 +0000 (05:25 +0000)
git-svn-id: https://svn.dealii.org/trunk@15314 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/numerics/data_component_interpretation.h [new file with mode: 0644]
deal.II/deal.II/include/numerics/data_out.h
deal.II/deal.II/include/numerics/data_postprocessor.h
deal.II/deal.II/source/numerics/data_out.cc
deal.II/deal.II/source/numerics/data_postprocessor.cc

diff --git a/deal.II/deal.II/include/numerics/data_component_interpretation.h b/deal.II/deal.II/include/numerics/data_component_interpretation.h
new file mode 100644 (file)
index 0000000..6e79cc0
--- /dev/null
@@ -0,0 +1,78 @@
+//---------------------------------------------------------------------------
+//    $Id: data_out.h 15313 2007-10-14 01:28:54Z bangerth $
+//    Version: $Name$
+//
+//    Copyright (C) 2007 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//---------------------------------------------------------------------------
+#ifndef __deal2__data_component_interpretation_h
+#define __deal2__data_component_interpretation_h
+
+
+
+#include <base/config.h>
+
+DEAL_II_NAMESPACE_OPEN
+
+/**
+ * A namespace solely for the declaration of the
+ * DataComponentInterpretation::DataComponentInterpretation enum.
+ */
+namespace DataComponentInterpretation
+{
+                                  /**
+                                   * The members of this enum are used to
+                                   * describe the logical interpretation of
+                                   * what the various components of a
+                                   * vector-valued data set mean. For
+                                   * example, if one has a finite element
+                                   * for the Stokes equations in 2d,
+                                   * representing components (u,v,p), one
+                                   * would like to indicate that the first
+                                   * two, u and v, represent a logical
+                                   * vector so that later on when we
+                                   * generate graphical output we can hand
+                                   * them off to a visualization program
+                                   * that will automatically know to render
+                                   * them as a vector field, rather than as
+                                   * two separate and independent scalar
+                                   * fields.
+                                   *
+                                   * By passing a set of enums of the
+                                   * current kind to the
+                                   * DataOut_DoFData::add_data_vector
+                                   * functions, this can be achieved.
+                                   *
+                                   * See the @ref step_22 "step-22" tutorial
+                                   * program for an example on how this
+                                   * information can be used in practice.
+                                   *
+                                   * @author Wolfgang Bangerth, 2007
+                                   */
+  enum DataComponentInterpretation
+  {
+                                        /**
+                                         * Indicates that a component of a
+                                         * data set corresponds to a scalar
+                                         * field independent of the others.
+                                         */
+       component_is_scalar,
+
+                                        /**
+                                         * Indicates that a component of a
+                                         * data set is part of a
+                                         * vector-valued quantity.
+                                         */
+       component_is_part_of_vector
+  };
+}
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
index 069d0ea65af5c1325851488ba48c8e1e226c2e78..fe22abb0728bf2540879df847e32dbc3c04de879 100644 (file)
@@ -22,6 +22,7 @@
 #include <dofs/dof_handler.h>
 #include <fe/mapping.h>
 #include <numerics/data_postprocessor.h>
+#include <numerics/data_component_interpretation.h>
 
 #include <boost/shared_ptr.hpp>
 
@@ -85,19 +86,28 @@ template <int> class FEValuesBase;
  * you can use the second add_data_vector() function which takes a @p
  * string instead of the <tt>vector<string></tt>.
  *
- * You should note that this class does not copy the vector given to it through
+ * The add_data_vector() functions have additional arguments (with default
+ * values) that can be used to specify certain transformations. In particular,
+ * it allows to attach DataPostprocessor arguments to compute derived
+ * information from a data vector at each quadrature point (for example, the
+ * Mach number in hypersonic flow can be computed from density and velocities;
+ * @ref step_29 "step-29" also shows an example); another piece of information
+ * specified through arguments with default values is how certain output
+ * components should be interpreted, i.e. whether each component of the data
+ * is logically an independent scalar field, or whether some of them together
+ * form logically a vector-field (see the
+ * DataComponentInterpretation::DataComponentInterpretation enum, and the @ref
+ * step_22 "step-22" tutorial program).
+ *
+ * It should be noted that this class does not copy the vector given to it through
  * the add_data_vector() functions, for memory consumption reasons. It only
  * stores a reference to it, so it is in your responsibility to make sure that
  * the data vectors exist long enough.
  *
  * After adding all data vectors, you need to call a function which generates
- * the patches for output from the stored data. This function is here called
- * build_patches(), but the naming is up to the derived class that actually
- * implements this.
- *
- * Finally, you write() the data in one format or other,
- * to a file and may want to clear this object as soon as possible to reduce
- * memory requirements.
+ * the patches for output from the stored data. Derived classes name this
+ * function build_patches(). Finally, you write() the data in one format or other,
+ * to a file.
  *
  * Please note, that in the example above, an object of type DataOut was
  * used, i.e. an object of a derived class. This is necessary since this
@@ -217,46 +227,6 @@ class DataOut_DoFData : public DataOutInterface<patch_dim,patch_space_dim>
                                            */
          type_automatic
     };
-
-                                    /**
-                                     * The members of this enum are used to
-                                     * describe the logical interpretation of
-                                     * what the various components of a
-                                     * vector-valued data set mean. For
-                                     * example, if one has a finite element
-                                     * for the Stokes equations in 2d,
-                                     * representing components (u,v,p), one
-                                     * would like to indicate that the first
-                                     * two, u and v, represent a logical
-                                     * vector so that later on when we
-                                     * generate graphical output we can hand
-                                     * them off to a visualization program
-                                     * that will automatically know to render
-                                     * them as a vector field, rather than as
-                                     * two separate and independent scalar
-                                     * fields.
-                                     *
-                                     * By passing a set of enums of the
-                                     * current kind to the
-                                     * DataOut_DoFData::add_data_vector
-                                     * functions, this can be achieved.
-                                     */
-    enum DataComponentInterpretation
-    {
-                                          /**
-                                           * Indicates that a component of a
-                                           * data set corresponds to a scalar
-                                           * field independent of the others.
-                                           */
-         component_is_scalar,
-
-                                          /**
-                                           * Indicates that a component of a
-                                           * data set is part of a
-                                           * vector-valued quantity.
-                                           */
-         component_is_part_of_vector
-    };
     
                                     /**
                                      * Constructor
@@ -340,6 +310,42 @@ class DataOut_DoFData : public DataOutInterface<patch_dim,patch_space_dim>
                                      * name instead of a vector of
                                      * names.
                                      *
+                                     * The data_component_interpretation
+                                     * argument contains information about
+                                     * how the individual components of
+                                     * output files that consist of more than
+                                     * one data set are to be interpreted.
+                                     *
+                                     * For example, if one has a finite
+                                     * element for the Stokes equations in
+                                     * 2d, representing components (u,v,p),
+                                     * one would like to indicate that the
+                                     * first two, u and v, represent a
+                                     * logical vector so that later on when
+                                     * we generate graphical output we can
+                                     * hand them off to a visualization
+                                     * program that will automatically know
+                                     * to render them as a vector field,
+                                     * rather than as two separate and
+                                     * independent scalar fields.
+                                     *
+                                     * The default value of this argument
+                                     * (i.e. an empty vector) corresponds is
+                                     * equivalent to a vector of values
+                                     * DataComponentInterpretation::component_is_scalar,
+                                     * indicating that all output components
+                                     * are independent scalar
+                                     * fields. However, if the given data
+                                     * vector represents logical vectors, you
+                                     * may pass a vector that contains values
+                                     * DataComponentInterpretation::component_is_part_of_vector. In
+                                     * the example above, one would pass in a
+                                     * vector with components
+                                     * (DataComponentInterpretation::component_is_part_of_vector,
+                                     * DataComponentInterpretation::component_is_part_of_vector,
+                                     * DataComponentInterpretation::component_is_scalar)
+                                     * for (u,v,p).
+                                     *
                                      * The names of a data vector
                                      * shall only contain characters
                                      * which are letters, underscore
@@ -358,23 +364,23 @@ class DataOut_DoFData : public DataOutInterface<patch_dim,patch_space_dim>
     template <class VECTOR>
     void add_data_vector (const VECTOR                   &data,
                          const std::vector<std::string> &names,
-                         const DataVectorType            type = type_automatic);
+                         const DataVectorType            type = type_automatic,
+                         const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation
+                         = std::vector<DataComponentInterpretation::DataComponentInterpretation>());
 
                                     /**
-                                     * This function is an
-                                     * abbreviation to the above one,
-                                     * intended for use with finite
-                                     * elements that are not composed
-                                     * of subelements. In this case,
-                                     * only one name per data vector
-                                     * needs to be given, which is
-                                     * what this function takes. It
-                                     * simply relays its arguments
-                                     * after a conversion of the
-                                     * @p name to a vector of
-                                     * strings, to the other
-                                     * add_data_vector() function
-                                     * above.
+                                     * This function is an abbreviation to
+                                     * the above one (see there for a
+                                     * discussion of the various arguments),
+                                     * intended for use with finite elements
+                                     * that are not composed of
+                                     * subelements. In this case, only one
+                                     * name per data vector needs to be
+                                     * given, which is what this function
+                                     * takes. It simply relays its arguments
+                                     * after a conversion of the @p name to a
+                                     * vector of strings, to the other
+                                     * add_data_vector() function above.
                                      *
                                      * If @p data is a vector with
                                      * multiple components this
@@ -393,7 +399,9 @@ class DataOut_DoFData : public DataOutInterface<patch_dim,patch_space_dim>
     template <class VECTOR>
     void add_data_vector (const VECTOR         &data,
                          const std::string    &name,
-                         const DataVectorType  type = type_automatic);
+                         const DataVectorType  type = type_automatic,
+                         const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation
+                         = std::vector<DataComponentInterpretation::DataComponentInterpretation>());
 
                                     /**
                                      * This function is an alternative to the
@@ -403,12 +411,16 @@ class DataOut_DoFData : public DataOutInterface<patch_dim,patch_space_dim>
                                      * a class derived from DataPostprocessor.
                                      *
                                      * The names for these derived quantities
-                                     * are implicitly given via the @p
-                                     * data_postprocessor argument. As only
-                                     * data of type @p type_dof_data can be
-                                     * transformed, this type is also known
-                                     * implicitly and does not have to be
-                                     * given.
+                                     * are provided by the @p
+                                     * data_postprocessor argument. Likewise,
+                                     * the data_component_interpretation
+                                     * argument of the other
+                                     * add_data_vector() functions is
+                                     * provided by the data_postprocessor
+                                     * argument. As only data of type @p
+                                     * type_dof_data can be transformed, this
+                                     * type is also known implicitly and does
+                                     * not have to be given.
                                      *
                                      * The actual type for the template
                                      * argument may be any vector type from
@@ -623,11 +635,14 @@ class DataOut_DoFData : public DataOutInterface<patch_dim,patch_space_dim>
                                         /**
                                          * Constructor. Give a list of names
                                          * for the individual components of
-                                         * the vector. This constructor
-                                         * assumes that no postprocessor is
-                                         * going to be used.
+                                         * the vector and their
+                                         * interpretation as scalar or vector
+                                         * data. This constructor assumes
+                                         * that no postprocessor is going to
+                                         * be used.
                                          */
-       DataEntryBase (const std::vector<std::string> &names);
+       DataEntryBase (const std::vector<std::string> &names,
+                      const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation);
 
                                         /**
                                          * Constructor when a data
@@ -751,7 +766,7 @@ class DataOut_DoFData : public DataOutInterface<patch_dim,patch_space_dim>
                                          * parts of a vector-field, or any of
                                          * the other supported kinds of data.
                                          */
-       const std::vector<DataComponentInterpretation>
+       const std::vector<DataComponentInterpretation::DataComponentInterpretation>
        data_component_interpretation;
        
                                         /**
@@ -790,12 +805,15 @@ class DataOut_DoFData : public DataOutInterface<patch_dim,patch_space_dim>
                                         /**
                                          * Constructor. Give a list of names
                                          * for the individual components of
-                                         * the vector. This constructor
-                                         * assumes that no postprocessor is
-                                         * going to be used.
+                                         * the vector and their
+                                         * interpretation as scalar or vector
+                                         * data. This constructor assumes
+                                         * that no postprocessor is going to
+                                         * be used.
                                          */
        DataEntry (const VectorType               *data,
-                  const std::vector<std::string> &names);
+                  const std::vector<std::string> &names,
+                  const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation);
 
                                         /**
                                          * Constructor when a data
@@ -965,8 +983,9 @@ class DataOut_DoFData : public DataOutInterface<patch_dim,patch_space_dim>
 
                                     /**
                                      * Overload of the respective
-                                     * DataOutBase::get_vector_data_ranges()
-                                     * function.
+                                     * DataOutInterface::get_vector_data_ranges()
+                                     * function. See there for a more
+                                     * extensive documentation.
                                      */
     virtual
     std::vector<boost::tuple<unsigned int, unsigned int, std::string> >
index 27639bb569d92c72c0c169d4ab14d4b1e578d14d..09251c52e8abce595f5766627322958c98d06c5e 100644 (file)
@@ -19,6 +19,7 @@
 #include <base/tensor.h>
 #include <lac/vector.h>
 #include <fe/fe_update_flags.h>
+#include <numerics/data_component_interpretation.h>
 
 #include <vector>
 #include <string>
@@ -135,7 +136,49 @@ class DataPostprocessor: public Subscriptor
                                      * the names of the computed quantities.
                                      */
     virtual std::vector<std::string> get_names () const=0;
-
+    
+                                    /**
+                                     * This functions returns
+                                     * information about how the
+                                     * individual components of
+                                     * output files that consist of
+                                     * more than one data set are to
+                                     * be interpreted.
+                                     *
+                                     * For example, if one has a finite
+                                     * element for the Stokes equations in
+                                     * 2d, representing components (u,v,p),
+                                     * one would like to indicate that the
+                                     * first two, u and v, represent a
+                                     * logical vector so that later on when
+                                     * we generate graphical output we can
+                                     * hand them off to a visualization
+                                     * program that will automatically know
+                                     * to render them as a vector field,
+                                     * rather than as two separate and
+                                     * independent scalar fields.
+                                     *
+                                     * The default implementation of this
+                                     * function returns a vector of values
+                                     * DataComponentInterpretation::component_is_scalar,
+                                     * indicating that all output components
+                                     * are independent scalar
+                                     * fields. However, if a derived class
+                                     * produces data that represents vectors,
+                                     * it may return a vector that contains
+                                     * values
+                                     * DataComponentInterpretation::component_is_part_of_vector. In
+                                     * the example above, one would return a
+                                     * vector with components
+                                     * (DataComponentInterpretation::component_is_part_of_vector,
+                                     * DataComponentInterpretation::component_is_part_of_vector,
+                                     * DataComponentInterpretation::component_is_scalar)
+                                     * for (u,v,p).
+                                     */
+    virtual
+    std::vector<DataComponentInterpretation::DataComponentInterpretation>
+    get_data_component_interpretation () const;
+    
                                     /**
                                      * Return, which data has to be provided to
                                      * compute the derived quantities. This has
index 503d5208beee1a083dc2a83d383279a35500c507..437dfca8334d10ae1d8d9c64c396c30f8038d84b 100644 (file)
@@ -37,11 +37,11 @@ DEAL_II_NAMESPACE_OPEN
 
 template <class DH, int patch_dim, int patch_space_dim>
 DataOut_DoFData<DH,patch_dim,patch_space_dim>::
-DataEntryBase::DataEntryBase (const std::vector<std::string> &names_in)
+DataEntryBase::DataEntryBase (const std::vector<std::string> &names_in,
+                             const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation)
                :
                names(names_in),
-               data_component_interpretation (names.size(),
-                                              component_is_scalar),
+               data_component_interpretation (data_component_interpretation),
                postprocessor(0),
                n_output_variables(names.size())
 {
@@ -69,8 +69,7 @@ DataOut_DoFData<DH,patch_dim,patch_space_dim>::
 DataEntryBase::DataEntryBase (const DataPostprocessor<DH::dimension> *data_postprocessor)
                :
                names(data_postprocessor->get_names()),
-               data_component_interpretation (names.size(),
-                                              component_is_scalar),
+               data_component_interpretation (data_postprocessor->get_data_component_interpretation()),
                postprocessor(data_postprocessor),
                n_output_variables(data_postprocessor->n_output_variables())
 {
@@ -114,9 +113,10 @@ template <typename VectorType>
 DataOut_DoFData<DH,patch_dim,patch_space_dim>::
 DataEntry<VectorType>::
 DataEntry (const VectorType                       *data,
-          const std::vector<std::string>         &names)
+          const std::vector<std::string>         &names,
+          const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation)
                :
-                DataEntryBase (names),
+                DataEntryBase (names, data_component_interpretation),
                vector (data)
 {}
 
@@ -291,7 +291,8 @@ void
 DataOut_DoFData<DH,patch_dim,patch_space_dim>::
 add_data_vector (const VECTOR                             &vec,
                 const std::string                        &name,
-                const DataVectorType                      type)
+                const DataVectorType                      type,
+                const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation)
 {
   Assert (dofs != 0, ExcNoDoFHandlerSelected ());
   const unsigned int n_components = dofs->get_fe().n_components ();
@@ -318,7 +319,7 @@ add_data_vector (const VECTOR                             &vec,
        }
     }
   
-  add_data_vector (vec, names, type);
+  add_data_vector (vec, names, type, data_component_interpretation);
 }
 
 
@@ -330,10 +331,20 @@ void
 DataOut_DoFData<DH,patch_dim,patch_space_dim>::
 add_data_vector (const VECTOR                             &vec,
                 const std::vector<std::string>           &names,
-                const DataVectorType                      type)
+                const DataVectorType                      type,
+                const std::vector<DataComponentInterpretation::DataComponentInterpretation> &data_component_interpretation_)
 {
   Assert (dofs != 0, ExcNoDoFHandlerSelected ());
 
+  const std::vector<DataComponentInterpretation::DataComponentInterpretation> &
+    data_component_interpretation
+    = (data_component_interpretation_.size() != 0
+       ?
+       data_component_interpretation_
+       :
+       std::vector<DataComponentInterpretation::DataComponentInterpretation>
+       (names.size(), DataComponentInterpretation::component_is_scalar));
+  
                                   // either cell data and one name,
                                   // or dof data and n_components names
   DataVectorType actual_type = type;
@@ -371,7 +382,8 @@ add_data_vector (const VECTOR                             &vec,
            Assert (false, ExcInternalError());
     }
 
-  DataEntryBase * new_entry = new DataEntry<VECTOR>(&vec, names);
+  DataEntryBase * new_entry = new DataEntry<VECTOR>(&vec, names,
+                                                   data_component_interpretation);
   if (actual_type == type_dof_data)
     dof_data.push_back (boost::shared_ptr<DataEntryBase>(new_entry));
   else
@@ -511,7 +523,7 @@ DataOut_DoFData<DH,patch_dim,patch_space_dim>::get_vector_data_ranges () const
                                       // the current function all we care
                                       // about is vector data
       if ((*d)->data_component_interpretation[i] ==
-         component_is_part_of_vector)
+         DataComponentInterpretation::component_is_part_of_vector)
        {
                                           // ensure that there is a
                                           // continuous number of next
@@ -524,7 +536,7 @@ DataOut_DoFData<DH,patch_dim,patch_space_dim>::get_vector_data_ranges () const
          for (unsigned int dd=1; dd<patch_space_dim; ++dd)
            Assert ((*d)->data_component_interpretation[i+dd]
                    ==
-                   component_is_part_of_vector,
+                   DataComponentInterpretation::component_is_part_of_vector,
                    ExcInvalidVectorDeclaration (i,
                                                 (*d)->names[i]));
 
@@ -1085,13 +1097,15 @@ template void \
 DataOut_DoFData<DH,D2,D3>:: \
 add_data_vector<V > (const V             &, \
                      const std::string   &, \
-                     const DataVectorType); \
+                     const DataVectorType,  \
+                    const std::vector<DataComponentInterpretation::DataComponentInterpretation> &); \
 \
 template void \
 DataOut_DoFData<DH,D2,D3>:: \
 add_data_vector<V > (const V                        &, \
                      const std::vector<std::string> &, \
-                     const DataVectorType); \
+                     const DataVectorType,  \
+                    const std::vector<DataComponentInterpretation::DataComponentInterpretation> &); \
 \
 template void \
 DataOut_DoFData<DH,D2,D3>:: \
index 61ea0bd6bb8b862c904bd2e7eda7e56efd36eb29..217b201479431e4c752a439351afced92f1392e9 100644 (file)
@@ -49,6 +49,22 @@ compute_derived_quantities_vector (const std::vector<Vector<double> >
   AssertThrow(false,ExcPureFunctionCalled());
 }
 
+
+
+template <int dim>
+std::vector<DataComponentInterpretation::DataComponentInterpretation>
+DataPostprocessor<dim>::get_data_component_interpretation () const
+{
+                                  // default implementation assumes that all
+                                  // components are independent scalars
+  return
+    std::vector<DataComponentInterpretation::DataComponentInterpretation>
+    (n_output_variables(),
+     DataComponentInterpretation::component_is_scalar);
+}
+
+
+
 // explicit instantiation
 template class DataPostprocessor<deal_II_dimension>;
 

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.