--- /dev/null
+//---------------------------------------------------------------------------
+// $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
#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>
* 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
*/
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
* 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
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
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
* 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
/**
* 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
* 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;
/**
/**
* 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
/**
* 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> >
#include <base/tensor.h>
#include <lac/vector.h>
#include <fe/fe_update_flags.h>
+#include <numerics/data_component_interpretation.h>
#include <vector>
#include <string>
* 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
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())
{
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())
{
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)
{}
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 ();
}
}
- add_data_vector (vec, names, type);
+ add_data_vector (vec, names, type, data_component_interpretation);
}
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;
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
// 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
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]));
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>:: \