* add_data_vector() for the
* method used).
*/
- enum DataVectorType {
- /// Data vector entries are associated to degrees of freedom
+ enum DataVectorType
+ {
+ /**
+ * Data vector entries are
+ * associated to degrees of freedom
+ */
type_dof_data,
- /// Data vector entries are one per grid cell
+
+ /**
+ * Data vector entries are one per
+ * grid cell
+ */
type_cell_data,
- /// Find out automatically
+
+ /**
+ * Find out automatically
+ */
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
* Exception
*/
DeclException0 (ExcIncompatiblePatchLists);
+
+ DeclException2 (ExcInvalidVectorDeclaration,
+ int, std::string,
+ << "When declaring that a number of components in a data\n"
+ << "set to be output logically form a vector instead of\n"
+ << "simply a set of scalar fields, you need to specify\n"
+ << "this for all relevant components. Furthermore,\n"
+ << "vectors must always consist of exactly <dim>\n"
+ << "components. However, the vector component at\n"
+ << "position " << arg1 << " with name <" << arg2
+ << " does not satisfy these conditions.");
protected:
/**
*/
const std::vector<std::string> names;
+ /**
+ * A vector that for each of the
+ * n_output_variables variables of
+ * the current data set indicates
+ * whether they are scalar fields,
+ * parts of a vector-field, or any of
+ * the other supported kinds of data.
+ */
+ const std::vector<DataComponentInterpretation>
+ data_component_interpretation;
+
/**
* Pointer to a DataPostprocessing
* object which shall be applied to
* classes.
*/
std::vector<Patch> patches;
-
+
/**
* Function by which the base
* class's functions get to know
* what patches they shall write
* to a file.
*/
- virtual const std::vector<Patch> &
- get_patches () const;
+ virtual
+ const std::vector<Patch> & get_patches () const;
/**
* Virtual function through
* obtained by the output functions
* of the base class.
*/
- virtual std::vector<std::string> get_dataset_names () const;
+ virtual
+ std::vector<std::string> get_dataset_names () const;
+ /**
+ * Overload of the respective
+ * DataOutBase::get_vector_data_ranges()
+ * function.
+ */
+ virtual
+ std::vector<boost::tuple<unsigned int, unsigned int, std::string> >
+ get_vector_data_ranges () const;
+
/**
* Make all template siblings
* friends. Needed for the
// names
Assert (get_dataset_names() == source.get_dataset_names(),
ExcIncompatibleDatasetNames());
- // make sure patches are compatible
+ // make sure patches are compatible. we'll
+ // assume that if the first respective
+ // patches are ok that all the other ones
+ // are ok as well
Assert (patches[0].n_subdivisions == source_patches[0].n_subdivisions,
ExcIncompatiblePatchLists());
Assert (patches[0].data.n_rows() == source_patches[0].data.n_rows(),
Assert (patches[0].data.n_cols() == source_patches[0].data.n_cols(),
ExcIncompatiblePatchLists());
+ // check equality of the vector data
+ // specifications
+ Assert (get_vector_data_ranges().size() ==
+ source.get_vector_data_ranges().size(),
+ ExcMessage ("Both sources need to declare the same components "
+ "as vectors."));
+ for (unsigned int i=0; i<get_vector_data_ranges().size(); ++i)
+ {
+ Assert (get_vector_data_ranges()[i].template get<0>() ==
+ source.get_vector_data_ranges()[i].template get<0>(),
+ ExcMessage ("Both sources need to declare the same components "
+ "as vectors."));
+ Assert (get_vector_data_ranges()[i].template get<1>() ==
+ source.get_vector_data_ranges()[i].template get<1>(),
+ ExcMessage ("Both sources need to declare the same components "
+ "as vectors."));
+ Assert (get_vector_data_ranges()[i].template get<2>() ==
+ source.get_vector_data_ranges()[i].template get<2>(),
+ ExcMessage ("Both sources need to declare the same components "
+ "as vectors."));
+ }
+
// merge patches. store old number
// of elements, since we need to
// adjust patch numbers, etc
const DataPostprocessor<DH::dimension> *data_postprocessor)
:
names(names_in),
+ data_component_interpretation (names.size(),
+ component_is_scalar),
postprocessor(data_postprocessor),
n_output_variables(names.size())
{
+//TODO this function should check that the number of names we have here equals the number of components to be obtained from the postprocessor
Assert(!data_postprocessor ||
data_postprocessor->n_output_variables()==names.size(),
- ExcDimensionMismatch(data_postprocessor->n_output_variables(),names.size()));
+ ExcDimensionMismatch(data_postprocessor->n_output_variables(),
+ names.size()));
+ Assert (!data_postprocessor ||
+ names.size() == data_component_interpretation.size(),
+ ExcDimensionMismatch(data_component_interpretation.size(),
+ names.size()));
}
{
Assert ((*d)->names.size() == 1, ExcInternalError());
names.push_back ((*d)->names[0]);
- };
+ }
return names;
}
+template <class DH,
+ int patch_dim, int patch_space_dim>
+std::vector<boost::tuple<unsigned int, unsigned int, std::string> >
+DataOut_DoFData<DH,patch_dim,patch_space_dim>::get_vector_data_ranges () const
+{
+ std::vector<boost::tuple<unsigned int, unsigned int, std::string> >
+ ranges;
+
+ // collect the ranges of dof
+ // and cell data
+ typedef
+ typename std::vector<boost::shared_ptr<DataEntryBase> >::const_iterator
+ data_iterator;
+
+ unsigned int output_component = 0;
+ for (data_iterator d=dof_data.begin();
+ d!=dof_data.end(); ++d)
+ for (unsigned int i=0; i<(*d)->n_output_variables;
+ ++i, ++output_component)
+ // see what kind of data we have
+ // here. note that for the purpose of
+ // the current function all we care
+ // about is vector data
+ if ((*d)->data_component_interpretation[i] ==
+ component_is_part_of_vector)
+ {
+ // ensure that there is a
+ // continuous number of next
+ // space_dim components that all
+ // deal with vectors
+ Assert (i+patch_space_dim <=
+ (*d)->n_output_variables,
+ ExcInvalidVectorDeclaration (i,
+ (*d)->names[i]));
+ for (unsigned int dd=1; dd<patch_space_dim; ++dd)
+ Assert ((*d)->data_component_interpretation[i+dd]
+ ==
+ component_is_part_of_vector,
+ ExcInvalidVectorDeclaration (i,
+ (*d)->names[i]));
+
+ // all seems alright, so figure out
+ // whether there is a common name
+ // to these components. if not,
+ // leave the name empty and let the
+ // output format writer decide what
+ // to do here
+ std::string name = (*d)->names[i];
+ for (unsigned int dd=1; dd<patch_space_dim; ++dd)
+ if (name != (*d)->names[i+dd])
+ {
+ name = "";
+ break;
+ }
+
+ // finally add a corresponding
+ // range
+ boost::tuple<unsigned int, unsigned int, std::string>
+ range (output_component,
+ output_component+patch_space_dim-1,
+ name);
+
+ ranges.push_back (range);
+
+ // increase the 'component' counter
+ // by the appropriate amount
+ output_component += patch_space_dim-1;
+ }
+
+ // note that we do not have to traverse the
+ // list of cell data here because cell data
+ // is one value per (logical) cell and
+ // therefore cannot be a vector
+
+ // as a final check, the 'component'
+ // counter should be at the total number of
+ // components added up now
+#ifdef DEBUG
+ unsigned int n_output_components = 0;
+ for (data_iterator d=dof_data.begin();
+ d!=dof_data.end(); ++d)
+ n_output_components = (*d)->n_output_variables;
+ Assert (output_component == n_output_components,
+ ExcInternalError());
+#endif
+
+ return ranges;
+}
+
+
+
template <class DH,
int patch_dim, int patch_space_dim>
const std::vector< dealii::DataOutBase::Patch<patch_dim, patch_space_dim> > &