Reference documentation for deal.II version Git 2e2d57a 2014-11-27 22:08:57 +0100
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Public Types | Public Member Functions | Private Member Functions | List of all members
DataOut< dim, DH > Class Template Reference

#include <deal.II/numerics/data_out.h>

Inheritance diagram for DataOut< dim, DH >:
[legend]

Public Types

enum  CurvedCellRegion { no_curved_cells, curved_boundary, curved_inner_cells }
 
typedef DataOut_DoFData< DH,
DH::dimension,
DH::space_dimension >
::cell_iterator 
cell_iterator
 
typedef DataOut_DoFData< DH,
DH::dimension,
DH::space_dimension >
::active_cell_iterator 
active_cell_iterator
 
- Public Types inherited from DataOut_DoFData< DH, DH::dimension, DH::space_dimension >
enum  DataVectorType
 
typedef Triangulation
< DH::dimension,
DH::space_dimension >
::cell_iterator 
cell_iterator
 
typedef Triangulation
< DH::dimension,
DH::space_dimension >
::active_cell_iterator 
active_cell_iterator
 

Public Member Functions

virtual void build_patches (const unsigned int n_subdivisions=0)
 
virtual void build_patches (const Mapping< DH::dimension, DH::space_dimension > &mapping, const unsigned int n_subdivisions=0, const CurvedCellRegion curved_region=curved_boundary)
 
virtual cell_iterator first_cell ()
 
virtual cell_iterator next_cell (const cell_iterator &cell)
 
 DeclException1 (ExcInvalidNumberOfSubdivisions, int,<< "The number of subdivisions per patch, "<< arg1<< ", is not valid.")
 
- Public Member Functions inherited from DataOut_DoFData< DH, DH::dimension, DH::space_dimension >
 DataOut_DoFData ()
 
virtual ~DataOut_DoFData ()
 
void attach_dof_handler (const DH &)
 
void attach_triangulation (const Triangulation< DH::dimension, DH::space_dimension > &)
 
void add_data_vector (const VECTOR &data, const std::vector< std::string > &names, const DataVectorType type=type_automatic, const std::vector< DataComponentInterpretation::DataComponentInterpretation > &data_component_interpretation=std::vector< DataComponentInterpretation::DataComponentInterpretation >())
 
void add_data_vector (const VECTOR &data, const std::string &name, const DataVectorType type=type_automatic, const std::vector< DataComponentInterpretation::DataComponentInterpretation > &data_component_interpretation=std::vector< DataComponentInterpretation::DataComponentInterpretation >())
 
void add_data_vector (const DH &dof_handler, const VECTOR &data, const std::vector< std::string > &names, const std::vector< DataComponentInterpretation::DataComponentInterpretation > &data_component_interpretation=std::vector< DataComponentInterpretation::DataComponentInterpretation >())
 
void add_data_vector (const DH &dof_handler, const VECTOR &data, const std::string &name, const std::vector< DataComponentInterpretation::DataComponentInterpretation > &data_component_interpretation=std::vector< DataComponentInterpretation::DataComponentInterpretation >())
 
void add_data_vector (const VECTOR &data, const DataPostprocessor< DH::space_dimension > &data_postprocessor)
 
void add_data_vector (const DH &dof_handler, const VECTOR &data, const DataPostprocessor< DH::space_dimension > &data_postprocessor)
 
void clear_data_vectors ()
 
void clear_input_data_references ()
 
void merge_patches (const DataOut_DoFData< DH2, patch_dim, patch_space_dim > &source, const Point< patch_space_dim > &shift=Point< patch_space_dim >())
 
virtual void clear ()
 
std::size_t memory_consumption () const
 
 DeclException0 (ExcNoTriangulationSelected)
 
 DeclException0 (ExcNoDoFHandlerSelected)
 
 DeclException0 (ExcDataPostprocessingIsNotPossibleForCellData)
 
 DeclException0 (ExcOldDataStillPresent)
 
 DeclException0 (ExcNoPatches)
 
 DeclException0 (ExcIncompatibleDatasetNames)
 
 DeclException0 (ExcIncompatiblePatchLists)
 
 DeclException3 (ExcInvalidVectorSize, int, int, int,<< "The vector has size "<< arg1<< " but the DoFHandler objects says there are "<< arg2<< " degrees of freedom and there are "<< arg3<< " active cells.")
 
 DeclException2 (ExcInvalidCharacter, std::string, size_t,<< "Please use only the characters [a-zA-Z0-9_<>()] for"<< std::endl<< "description strings since some graphics formats will only accept these."<< std::endl<< "The string you gave was <"<< arg1<< ">, the invalid character is <"<< arg1[arg2]<< ">."<< std::endl)
 
 DeclException2 (ExcInvalidNumberOfNames, int, int,<< "You have to give one name per component in your "<< "data vector. The number you gave was "<< arg1<< ", but the number of components is "<< arg2)
 
 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.")
 
- Public Member Functions inherited from DataOutInterface< patch_dim, patch_space_dim >
 DataOutInterface ()
 
virtual ~DataOutInterface ()
 
void write_dx (std::ostream &out) const
 
void write_eps (std::ostream &out) const
 
void write_gmv (std::ostream &out) const
 
void write_gnuplot (std::ostream &out) const
 
void write_povray (std::ostream &out) const
 
void write_tecplot (std::ostream &out) const
 
void write_tecplot_binary (std::ostream &out) const
 
void write_ucd (std::ostream &out) const
 
void write_vtk (std::ostream &out) const
 
void write_vtu (std::ostream &out) const
 
void write_vtu_in_parallel (const char *filename, MPI_Comm comm) const
 
void write_pvtu_record (std::ostream &out, const std::vector< std::string > &piece_names) const
 
void write_pvd_record (std::ostream &out, const std::vector< std::pair< double, std::string > > &times_and_names) const
 
void write_visit_record (std::ostream &out, const std::vector< std::string > &piece_names) const
 
void write_visit_record (std::ostream &out, const std::vector< std::vector< std::string > > &piece_names) const
 
void write_svg (std::ostream &out) const
 
void write_deal_II_intermediate (std::ostream &out) const
 
XDMFEntry create_xdmf_entry (const std::string &h5_filename, const double cur_time, MPI_Comm comm) const DEAL_II_DEPRECATED
 
XDMFEntry create_xdmf_entry (const DataOutBase::DataOutFilter &data_filter, const std::string &h5_filename, const double cur_time, MPI_Comm comm) const
 
XDMFEntry create_xdmf_entry (const DataOutBase::DataOutFilter &data_filter, const std::string &h5_mesh_filename, const std::string &h5_solution_filename, const double cur_time, MPI_Comm comm) const
 
void write_xdmf_file (const std::vector< XDMFEntry > &entries, const std::string &filename, MPI_Comm comm) const
 
void write_hdf5_parallel (const std::string &filename, MPI_Comm comm) const DEAL_II_DEPRECATED
 
void write_hdf5_parallel (const DataOutBase::DataOutFilter &data_filter, const std::string &filename, MPI_Comm comm) const
 
void write_hdf5_parallel (const DataOutBase::DataOutFilter &data_filter, const bool write_mesh_file, const std::string &mesh_filename, const std::string &solution_filename, MPI_Comm comm) const
 
void write_filtered_data (DataOutBase::DataOutFilter &filtered_data) const
 
void write (std::ostream &out, const DataOutBase::OutputFormat output_format=DataOutBase::default_format) const
 
void set_default_format (const DataOutBase::OutputFormat default_format)
 
void set_flags (const DataOutBase::DXFlags &dx_flags)
 
void set_flags (const DataOutBase::UcdFlags &ucd_flags)
 
void set_flags (const DataOutBase::GnuplotFlags &gnuplot_flags)
 
void set_flags (const DataOutBase::PovrayFlags &povray_flags)
 
void set_flags (const DataOutBase::EpsFlags &eps_flags)
 
void set_flags (const DataOutBase::GmvFlags &gmv_flags)
 
void set_flags (const DataOutBase::TecplotFlags &tecplot_flags)
 
void set_flags (const DataOutBase::VtkFlags &vtk_flags)
 
void set_flags (const DataOutBase::SvgFlags &svg_flags)
 
void set_flags (const DataOutBase::Deal_II_IntermediateFlags &deal_II_intermediate_flags)
 
std::string default_suffix (const DataOutBase::OutputFormat output_format=DataOutBase::default_format) const
 
void parse_parameters (ParameterHandler &prm)
 
std::size_t memory_consumption () const
 

Private Member Functions

cell_iterator first_locally_owned_cell ()
 
cell_iterator next_locally_owned_cell (const cell_iterator &cell)
 
void build_one_patch (const std::pair< cell_iterator, unsigned int > *cell_and_index, internal::DataOut::ParallelData< DH::dimension, DH::space_dimension > &data,::DataOutBase::Patch< DH::dimension, DH::space_dimension > &patch, const CurvedCellRegion curved_cell_region, std::vector<::DataOutBase::Patch< DH::dimension, DH::space_dimension > > &patches)
 

Additional Inherited Members

- Static Public Member Functions inherited from DataOutInterface< patch_dim, patch_space_dim >
static void declare_parameters (ParameterHandler &prm)
 
- Protected Types inherited from DataOut_DoFData< DH, DH::dimension, DH::space_dimension >
typedef ::DataOutBase::Patch
< patch_dim, patch_space_dim > 
Patch
 
- Protected Member Functions inherited from DataOut_DoFData< DH, DH::dimension, DH::space_dimension >
virtual const std::vector
< Patch > & 
get_patches () const
 
virtual std::vector< std::string > get_dataset_names () const
 
std::vector
< std_cxx11::shared_ptr
<::hp::FECollection
< DH::dimension,
DH::space_dimension > > > 
get_finite_elements () const
 
virtual std::vector
< std_cxx11::tuple< unsigned
int, unsigned int, std::string > > 
get_vector_data_ranges () const
 
- Protected Attributes inherited from DataOut_DoFData< DH, DH::dimension, DH::space_dimension >
SmartPointer< const
Triangulation< DH::dimension,
DH::space_dimension > > 
triangulation
 
SmartPointer< const DH > dofs
 
std::vector
< std_cxx11::shared_ptr
< internal::DataOut::DataEntryBase
< DH > > > 
dof_data
 
std::vector
< std_cxx11::shared_ptr
< internal::DataOut::DataEntryBase
< DH > > > 
cell_data
 
std::vector< Patchpatches
 
- Protected Attributes inherited from DataOutInterface< patch_dim, patch_space_dim >
unsigned int default_subdivisions
 

Detailed Description

template<int dim, class DH = DoFHandler<dim>>
class DataOut< dim, DH >

This class is the main class to provide output of data described by finite element fields defined on a collection of cells.

This class is an actual implementation of the functionality proposed by the DataOut_DoFData class. It offers a function build_patches() that generates the patches to be written in some graphics format from the data stored in the base class. Most of the interface and an example of its use is described in the documentation of the base class.

The only thing this class offers is the function build_patches() which loops over all cells of the triangulation stored by the attach_dof_handler() function of the base class (with the exception of cells of parallel::distributed::Triangulation objects that are not owned by the current processor) and converts the data on these to actual patches which are the objects that are later output by the functions of the base classes. You can give a parameter to the function which determines how many subdivisions in each coordinate direction are to be performed, i.e. of how many subcells each patch shall consist. Default is one, but you may want to choose a higher number for higher order elements, for example two for quadratic elements, three for cubic elements three, and so on. The purpose of this parameter is because most graphics programs do not allow to specify higher order polynomial functions in the file formats: only data at vertices can be plotted and is then shown as a bilinear interpolation within the interior of cells. This may be insufficient if you have higher order finite elements, and the only way to achieve better output is to subdivide each cell of the mesh into several cells for graphical output. Of course, what you get to see is still a bilinear interpolation on each cell of the output (where these cells are not subdivisions of the cells of the triangulation in use) due to the same limitations in output formats, but at least a bilinear interpolation of a higher order polynomial on a finer mesh.

Note that after having called build_patches() once, you can call one or more of the write() functions of DataOutInterface. You can therefore output the same data in more than one format without having to rebuild the patches.

User interface information

The base classes of this class, DataOutBase, DataOutInterface and DataOut_DoFData() offer several interfaces of their own. Refer to the DataOutBase class's documentation for a discussion of the different output formats presently supported, DataOutInterface for ways of selecting which format to use upon output at run-time and without the need to adapt your program when new formats become available, as well as for flags to determine aspects of output. The DataOut_DoFData() class's documentation has an example of using nodal data to generate output.

Extensions

By default, this class produces patches for all active cells. Sometimes, this is not what you want, maybe because they are simply too many (and too small to be seen individually) or because you only want to see a certain region of the domain (for example in parallel programs such as the step-18 example program), or for some other reason.

For this, internally build_patches() does not generate the sequence of cells to be converted into patches itself, but relies on the two functions first_cell() and next_cell(). By default, they return the first active cell, and the next active cell, respectively. Since they are virtual functions, you can write your own class derived from DataOut in which you overload these two functions to select other cells for output. This may, for example, include only cells that are in parts of a domain (e.g., if you don't care about the solution elsewhere, think for example a buffer region in which you attenuate outgoing waves in the Perfectly Matched Layer method) or if you don't want output to be generated at all levels of an adaptively refined mesh because this creates too much data (in this case, the set of cells returned by your implementations of first_cell() and next_cell() will include non-active cells, and DataOut::build_patches() will simply take interpolated values of the solution instead of the exact values on these cells children for output). Once you derive your own class, you would just create an object of this type instead of an object of type DataOut, and everything else will remain the same.

The two functions are not constant, so you may store information within your derived class about the last accessed cell. This is useful if the information of the last cell which was accessed is not sufficient to determine the next one.

There is one caveat, however: if you have cell data (in contrast to nodal, or dof, data) such as error indicators, then you must make sure that first_cell() and next_cell() only walk over active cells, since cell data cannot be interpolated to a coarser cell. If you do have cell data and use this pair of functions and they return a non-active cell, then an exception will be thrown.

Precondition
This class only makes sense if the first template argument, dim equals the dimension of the DoFHandler type given as the second template argument, i.e., if dim == DH::dimension. This redundancy is a historical relic from the time where the library had only a single DoFHandler class and this class consequently only a single template argument.
Author
Wolfgang Bangerth, 1999

Definition at line 161 of file data_out.h.

Member Typedef Documentation

template<int dim, class DH = DoFHandler<dim>>
typedef DataOut_DoFData<DH,DH::dimension,DH::space_dimension>::cell_iterator DataOut< dim, DH >::cell_iterator

Typedef to the iterator type of the dof handler class under consideration.

Definition at line 168 of file data_out.h.

Member Enumeration Documentation

template<int dim, class DH = DoFHandler<dim>>
enum DataOut::CurvedCellRegion

Enumeration describing the region of the domain in which curved cells shall be created.

Definition at line 175 of file data_out.h.

Member Function Documentation

template<int dim, class DH = DoFHandler<dim>>
virtual void DataOut< dim, DH >::build_patches ( const unsigned int  n_subdivisions = 0)
virtual

This is the central function of this class since it builds the list of patches to be written by the low-level functions of the base class. See the general documentation of this class for further information.

The default value 0 of n_subdivisions indicates that the value stored in DataOutInterface::default_subdivisions is to be used.

template<int dim, class DH = DoFHandler<dim>>
virtual void DataOut< dim, DH >::build_patches ( const Mapping< DH::dimension, DH::space_dimension > &  mapping,
const unsigned int  n_subdivisions = 0,
const CurvedCellRegion  curved_region = curved_boundary 
)
virtual

Same as above, except that the additional first parameter defines a mapping that is to be used in the generation of output. If n_subdivisions>1, the points interior of subdivided patches which originate from cells at the boundary of the domain can be computed using the mapping, i.e. a higher order mapping leads to a representation of a curved boundary by using more subdivisions. Some mappings like MappingQEulerian result in curved cells in the interior of the domain. However, there is nor easy way to get this information from the Mapping. Thus the last argument curved_region take one of three values resulting in no curved cells at all, curved cells at the boundary (default) or curved cells in the whole domain.

Even for non-curved cells the mapping argument can be used for the Eulerian mappings (see class MappingQ1Eulerian) where a mapping is used not only to determine the position of points interior to a cell, but also of the vertices. It offers an opportunity to watch the solution on a deformed triangulation on which the computation was actually carried out, even if the mesh is internally stored in its undeformed configuration and the deformation is only tracked by an additional vector that holds the deformation of each vertex.

Todo:
The mapping argument should be replaced by a hp::MappingCollection in case of a hp::DoFHandler.
template<int dim, class DH = DoFHandler<dim>>
virtual cell_iterator DataOut< dim, DH >::first_cell ( )
virtual

Return the first cell which we want output for. The default implementation returns the first active cell, but you might want to return other cells in a derived class.

template<int dim, class DH = DoFHandler<dim>>
virtual cell_iterator DataOut< dim, DH >::next_cell ( const cell_iterator cell)
virtual

Return the next cell after cell which we want output for. If there are no more cells, dofs->end() shall be returned.

The default implementation returns the next active cell, but you might want to return other cells in a derived class. Note that the default implementation assumes that the given cell is active, which is guaranteed as long as first_cell() is also used from the default implementation. Overloading only one of the two functions might not be a good idea.

template<int dim, class DH = DoFHandler<dim>>
DataOut< dim, DH >::DeclException1 ( ExcInvalidNumberOfSubdivisions  ,
int  ,
<< "The number of subdivisions per  patch,
"<< arg1<< "  ,
is not valid."   
)

Exception

template<int dim, class DH = DoFHandler<dim>>
cell_iterator DataOut< dim, DH >::first_locally_owned_cell ( )
private

Return the first cell produced by the first_cell()/next_cell() function pair that is locally owned. If this object operates on a non-distributed triangulation, the result equals what first_cell() returns.

template<int dim, class DH = DoFHandler<dim>>
cell_iterator DataOut< dim, DH >::next_locally_owned_cell ( const cell_iterator cell)
private

Return the next cell produced by the next_cell() function that is locally owned. If this object operates on a non-distributed triangulation, the result equals what first_cell() returns.

template<int dim, class DH = DoFHandler<dim>>
void DataOut< dim, DH >::build_one_patch ( const std::pair< cell_iterator, unsigned int > *  cell_and_index,
internal::DataOut< dim, DH >::ParallelData< DH::dimension, DH::space_dimension > &  data,
::DataOutBase::Patch< DH::dimension, DH::space_dimension > &  patch,
const CurvedCellRegion  curved_cell_region,
std::vector<::DataOutBase::Patch< DH::dimension, DH::space_dimension > > &  patches 
)
private

Build one patch. This function is called in a WorkStream context.

The result is written into the patch variable.


The documentation for this class was generated from the following file: