![]() |
Reference documentation for deal.II version 8.1.0
|
#include <data_out.h>
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 |
![]() | |
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 |
![]() | |
enum | OutputFormat |
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.") | |
![]() | |
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.") | |
![]() | |
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 > > ×_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 DataOutFilter &data_filter, const std::string &h5_filename, const double cur_time, MPI_Comm comm) const |
XDMFEntry | create_xdmf_entry (const 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 DataOutFilter &data_filter, const std::string &filename, MPI_Comm comm) const |
void | write_hdf5_parallel (const 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 (DataOutFilter &filtered_data) const |
void | write (std::ostream &out, const OutputFormat output_format=default_format) const |
void | set_default_format (const OutputFormat default_format) |
void | set_flags (const DXFlags &dx_flags) |
void | set_flags (const UcdFlags &ucd_flags) |
void | set_flags (const GnuplotFlags &gnuplot_flags) |
void | set_flags (const PovrayFlags &povray_flags) |
void | set_flags (const EpsFlags &eps_flags) |
void | set_flags (const GmvFlags &gmv_flags) |
void | set_flags (const TecplotFlags &tecplot_flags) |
void | set_flags (const VtkFlags &vtk_flags) |
void | set_flags (const SvgFlags &svg_flags) |
void | set_flags (const Deal_II_IntermediateFlags &deal_II_intermediate_flags) |
std::string | default_suffix (const OutputFormat output_format=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 void | declare_parameters (ParameterHandler &prm) |
![]() | |
typedef ::DataOutBase::Patch < patch_dim, patch_space_dim > | Patch |
![]() | |
virtual const std::vector < Patch > & | get_patches () const |
virtual std::vector< std::string > | get_dataset_names () const |
std::vector < std_cxx1x::shared_ptr <::hp::FECollection < DH::dimension, DH::space_dimension > > > | get_finite_elements () const |
virtual std::vector < std_cxx1x::tuple< unsigned int, unsigned int, std::string > > | get_vector_data_ranges () const |
![]() | |
SmartPointer< const Triangulation< DH::dimension, DH::space_dimension > > | triangulation |
SmartPointer< const DH > | dofs |
std::vector < std_cxx1x::shared_ptr < internal::DataOut::DataEntryBase < DH > > > | dof_data |
std::vector < std_cxx1x::shared_ptr < internal::DataOut::DataEntryBase < DH > > > | cell_data |
std::vector< Patch > | patches |
![]() | |
unsigned int | default_subdivisions |
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.
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.
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 may overload them to select other cells for output. If cells are not active, interpolated values are taken for output instead of the exact values on active cells.
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.
Definition at line 144 of file data_out.h.
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 151 of file data_out.h.
enum DataOut::CurvedCellRegion |
Enumeration describing the region of the domain in which curved cells shall be created.
Definition at line 158 of file data_out.h.
|
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.
|
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.
mapping
argument should be replaced by a hp::MappingCollection in case of a hp::DoFHandler.
|
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.
|
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.
DataOut< dim, DH >::DeclException1 | ( | ExcInvalidNumberOfSubdivisions | , |
int | , | ||
<< "The number of subdivisions per | patch, | ||
"<< arg1<< " | , | ||
is not valid." | |||
) |
Exception
|
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.
|
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.
|
private |
Build one patch. This function is called in a WorkStream context.
The result is written into the patch variable.