Reference documentation for deal.II version GIT relicensing-136-gb80d0be4af 2024-03-18 08:20:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
Classes | Enumerations | Functions
DataOutBase Namespace Reference

Classes

class  DataOutFilter
 
struct  DataOutFilterFlags
 
struct  Deal_II_IntermediateFlags
 
struct  DXFlags
 
struct  EpsFlags
 
struct  GmvFlags
 
struct  GnuplotFlags
 
struct  Hdf5Flags
 
struct  OutputFlagsBase
 
struct  Patch
 
struct  Patch< 0, spacedim >
 
struct  PovrayFlags
 
struct  SvgFlags
 
struct  TecplotFlags
 
struct  UcdFlags
 
struct  VtkFlags
 

Enumerations

enum class  CompressionLevel {
  no_compression , best_speed , best_compression , default_compression ,
  plain_text
}
 
enum  OutputFormat {
  default_format , none , dx , ucd ,
  gnuplot , povray , eps , gmv ,
  tecplot , vtk , vtu , svg ,
  deal_II_intermediate , hdf5
}
 

Functions

template<int dim, int spacedim>
void write_dx (const std::vector< Patch< dim, spacedim > > &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &nonscalar_data_ranges, const DXFlags &flags, std::ostream &out)
 
template<int spacedim>
void write_eps (const std::vector< Patch< 2, spacedim > > &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &nonscalar_data_ranges, const EpsFlags &flags, std::ostream &out)
 
template<int dim, int spacedim>
void write_eps (const std::vector< Patch< dim, spacedim > > &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &nonscalar_data_ranges, const EpsFlags &flags, std::ostream &out)
 
template<int dim, int spacedim>
void write_gmv (const std::vector< Patch< dim, spacedim > > &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &nonscalar_data_ranges, const GmvFlags &flags, std::ostream &out)
 
template<int dim, int spacedim>
void write_gnuplot (const std::vector< Patch< dim, spacedim > > &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &nonscalar_data_ranges, const GnuplotFlags &flags, std::ostream &out)
 
template<int dim, int spacedim>
void write_povray (const std::vector< Patch< dim, spacedim > > &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &nonscalar_data_ranges, const PovrayFlags &flags, std::ostream &out)
 
template<int dim, int spacedim>
void write_tecplot (const std::vector< Patch< dim, spacedim > > &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &nonscalar_data_ranges, const TecplotFlags &flags, std::ostream &out)
 
template<int dim, int spacedim>
void write_ucd (const std::vector< Patch< dim, spacedim > > &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &nonscalar_data_ranges, const UcdFlags &flags, std::ostream &out)
 
template<int dim, int spacedim>
void write_vtk (const std::vector< Patch< dim, spacedim > > &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &nonscalar_data_ranges, const VtkFlags &flags, std::ostream &out)
 
template<int dim, int spacedim>
void write_vtu (const std::vector< Patch< dim, spacedim > > &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &nonscalar_data_ranges, const VtkFlags &flags, std::ostream &out)
 
void write_vtu_header (std::ostream &out, const VtkFlags &flags)
 
void write_vtu_footer (std::ostream &out)
 
template<int dim, int spacedim>
void write_vtu_main (const std::vector< Patch< dim, spacedim > > &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &nonscalar_data_ranges, const VtkFlags &flags, std::ostream &out)
 
void write_pvtu_record (std::ostream &out, const std::vector< std::string > &piece_names, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &nonscalar_data_ranges, const VtkFlags &flags)
 
void write_pvd_record (std::ostream &out, const std::vector< std::pair< double, std::string > > &times_and_names)
 
void write_visit_record (std::ostream &out, const std::vector< std::string > &piece_names)
 
void write_visit_record (std::ostream &out, const std::vector< std::vector< std::string > > &piece_names)
 
void write_visit_record (std::ostream &out, const std::vector< std::pair< double, std::vector< std::string > > > &times_and_piece_names)
 
template<int spacedim>
void write_svg (const std::vector< Patch< 2, spacedim > > &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &nonscalar_data_ranges, const SvgFlags &flags, std::ostream &out)
 
template<int dim, int spacedim>
void write_deal_II_intermediate (const std::vector< Patch< dim, spacedim > > &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &nonscalar_data_ranges, const Deal_II_IntermediateFlags &flags, std::ostream &out)
 
template<int dim, int spacedim>
void write_deal_II_intermediate_in_parallel (const std::vector< Patch< dim, spacedim > > &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &nonscalar_data_ranges, const Deal_II_IntermediateFlags &flags, const std::string &filename, const MPI_Comm comm, const CompressionLevel compression)
 
template<int dim, int spacedim>
void write_hdf5_parallel (const std::vector< Patch< dim, spacedim > > &patches, const DataOutFilter &data_filter, const DataOutBase::Hdf5Flags &flags, const std::string &filename, const MPI_Comm comm)
 
template<int dim, int spacedim>
void write_hdf5_parallel (const std::vector< Patch< dim, spacedim > > &patches, const DataOutFilter &data_filter, const DataOutBase::Hdf5Flags &flags, const bool write_mesh_file, const std::string &mesh_filename, const std::string &solution_filename, const MPI_Comm comm)
 
template<int dim, int spacedim>
void write_filtered_data (const std::vector< Patch< dim, spacedim > > &patches, const std::vector< std::string > &data_names, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &nonscalar_data_ranges, DataOutFilter &filtered_data)
 
std::pair< unsigned int, unsigned intdetermine_intermediate_format_dimensions (std::istream &input)
 
OutputFormat parse_output_format (const std::string &format_name)
 
std::string get_output_format_names ()
 
std::string default_suffix (const OutputFormat output_format)
 
static ::ExceptionBaseExcInvalidDatasetSize (int arg1, int arg2)
 
static ::ExceptionBaseExcNoPatches ()
 
static ::ExceptionBaseExcTecplotAPIError ()
 
static ::ExceptionBaseExcErrorOpeningTecplotFile (char *arg1)
 
template<int dim, int spacedim>
std::ostream & operator<< (std::ostream &out, const Patch< dim, spacedim > &patch)
 
template<int dim, int spacedim>
std::istream & operator>> (std::istream &in, Patch< dim, spacedim > &patch)
 
template<int dim, int spacedim>
std::vector< Point< spacedim > > get_node_positions (const std::vector< Patch< dim, spacedim > > &patches)
 
template<int dim, int spacedim, typename StreamType >
void write_nodes (const std::vector< Patch< dim, spacedim > > &patches, StreamType &out)
 
template<int dim, int spacedim, typename StreamType >
void write_cells (const std::vector< Patch< dim, spacedim > > &patches, StreamType &out)
 
template<int dim, int spacedim, typename StreamType >
void write_high_order_cells (const std::vector< Patch< dim, spacedim > > &patches, StreamType &out, const bool legacy_format)
 
template<int dim, int spacedim, typename StreamType >
void write_data (const std::vector< Patch< dim, spacedim > > &patches, unsigned int n_data_sets, const bool double_precision, StreamType &out)
 
template<int dim, int spacedim>
void write_svg (const std::vector< Patch< dim, spacedim > > &, const std::vector< std::string > &, const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &, const SvgFlags &, std::ostream &)
 

Detailed Description

This is a base class for output of data on meshes of very general form. Output data is expected as a set of patches and written to the output stream in the format expected by the visualization tool. For a list of output formats, check the enumeration OutputFormat. For each format listed there, this class contains a function write_format, writing the output. Refer to the documentation of those functions for details on a certain format.

Structure of the output data

Data is not written with the deal.II mesh structure. Instead, it relies on a set of patches created by a derived class (for example the DataOut, DataOutStack, DataOutFaces, DataOutRotation, or MatrixOut classes). Each Patch describes a single logical cell of a mesh, possibly subdivided a number of times to represent higher order polynomials defined on this cell. To this end, a patch consists of a dim-dimensional regular grid with the same number of grid points in each direction. In the simplest case it may consist of the corner points of a single mesh cell. For each point of this local grid, the Patch contains an arbitrary number of data values, though the number of data sets must be the same for each point on each patch.

By offering this interface to the different output formats, it is simple to extend this class to new formats without depending on such things as actual triangulations and handling of data vectors. These things shall be provided by derived class which have a user callable interface then.

Inside each patch, the data is organized in the usual lexicographical order, x running fastest, then y and z. Nodes are stored in this order and cells as well. Each cell in 3d is stored such that the front face is in the xz-plane. In order to enhance intelligibility of this concept, the following two sections are kept from a previous version of this documentation.

Patches

Grids can be thought of as a collection of cells; if you want to write out data on such a grid, you can do so by writing them one cell at a time. The functions in this class therefore take a list of objects describing the data on one cell each. This data for each cell usually consists of a list of vertices for this cell, and a list of data values (for example solution data, error information, etc) at each of these vertices.

In some cases, this interface to a cell is too restricted, however. For example, you may have higher order elements and printing the values at the vertices only is not enough. For this reason, we not only provide writing the data on the vertices only, but the data is organizes as a tensor product grid on each cell. The parameter n_subdivisions, which is given for each patch separately, denotes how often the cell is to be divided for output; for example, n_subdivisions==1 yields no subdivision of the cell, n_subdivisions==2 will produce a grid of 3 times 3 points in two spatial dimensions and 3 times 3 times 3 points in three dimensions, n_subdivisions==3 will yield 4 times 4 (times 4) points, etc. The actual location of these points on the patch will be computed by a multilinear transformation from the vertices given for this patch. For cells at the boundary, a mapping might be used to calculate the position of the inner points. In that case the coordinates are stored inside the Patch, as they cannot be easily recovered otherwise.

Given these comments, the actual data to be printed on this patch of points consists of several data sets each of which has a value at each of the patch points. For example with n_subdivisions==2 in two space dimensions, each data set has to provide nine values, and since the patch is to be printed as a tensor product (or its transformation to the real space cell), its values are to be ordered like (x0,y0) (x0,y1) (x0,y2) (x1,y0) (x1,y1) (x1,y2) (x2,y0) (x2,y1) (x2,y2), i.e. the z-coordinate runs fastest, then the y-coordinate, then x (if there are that many space directions).

Generalized patches

In general, the patches as explained above might be too restricted. For example, one might want to draw only the outer faces of a domain in a three-dimensional computation, if one is not interested in what happens inside. Then, the objects that should be drawn are two-dimensional in a three-dimensional world. The Patch class and associated output functions handle these cases. The Patch class therefore takes two template parameters, the first, named dim denoting the dimension of the object (in the above example, this would be two), while the second, named spacedim, denotes the dimension of the embedding space (this would be three). The corner points of a patch have the dimension of the space, while their number is determined by the dimension of the patch. By default, the second template parameter has the same value as the first, which would correspond to outputting a cell, rather than a face or something else.

DataOutBaseInterface

The members of this namespace are not usually called from user code directly. Rather, classes that use the functions declared here are typically derived from DataOutInterface.

The interface of this class basically consists of the declaration of a data type describing a patch and a bunch of functions taking a list of patches and writing them in one format or other to the stream. It is in the responsibility of the derived classes to provide this list of patches. In addition to the list of patches, a name for each data set may be given.

Querying interface

This class also provides a few functions (parse_output_format(), get_output_format_names(), default_suffix()) that can be used to query which output formats this class supports. The provide a list of names for all the formats we can output, parse a string and return an enum indicating each format, and provide a way to convert a value of this enum into the usual suffix used for files of that name. Using these functions, one can entirely free applications from knowledge which formats the library presently allows to output; several of the example programs show how to do this.

Output parameters

All functions take a parameter which is a structure of type XFlags, where X is the name of the output format. To find out what flags are presently supported, read the documentation of the different structures.

Note that usually the output formats used for scientific visualization programs have no or very few parameters (apart from some compatibility flags) because there the actual appearance of output is determined using the visualization program and the files produced by this class store more or less only raw data.

The direct output formats, like Postscript or Povray need to be given a lot more parameters, though, since there the output file has to contain all details of the viewpoint, light source, etc.

Writing backends

An abstraction layer has been introduced to facilitate coding backends for additional visualization tools. It is applicable for data formats separating the information into a field of vertices, a field of connection information for the grid cells and data fields.

For each of these fields, output functions are implemented, namely write_nodes(), write_cells() and write_data(). In order to use these functions, a format specific output stream must be written, following the examples of DXStream, GmvStream, VtkStream and so on, implemented in the .cc file.

In this framework, the implementation of a new output format is reduced to writing the section headers and the new output stream class for writing a single mesh object.

Enumeration Type Documentation

◆ CompressionLevel

enum class DataOutBase::CompressionLevel
strong

An enum for different levels of compression used in several places to determine zlib compression levels for binary output. At some places, it is possible to output the data also as plain text as an alternative, which is convenient for debugging. We use this flag to indicate such an output as well.

Enumerator
no_compression 

Do not use any compression.

best_speed 

Use the fastest available compression algorithm.

best_compression 

Use the algorithm which results in the smallest compressed files.

default_compression 

Use the default compression algorithm. This is a compromise between speed and file size.

plain_text 

Output as plain text (ASCII) if available.

Definition at line 206 of file data_out_base.h.

◆ OutputFormat

Provide a data type specifying the presently supported output formats.

Enumerator
default_format 

Use the format already stored in the object.

none 

Do not write any output.

dx 

Output for OpenDX.

ucd 

Output in the UCD format for AVS.

gnuplot 

Output for the Gnuplot tool.

povray 

Output for the Povray raytracer.

eps 

Output in encapsulated PostScript.

gmv 

Output for GMV.

tecplot 

Output for Tecplot in text format.

vtk 

Output in VTK format.

vtu 

Output in VTK format.

svg 

Output in SVG format.

deal_II_intermediate 

Output in deal.II intermediate format.

hdf5 

Output in HDF5 format.

Definition at line 1612 of file data_out_base.h.

Function Documentation

◆ write_dx()

template<int dim, int spacedim>
void DataOutBase::write_dx ( const std::vector< Patch< dim, spacedim > > &  patches,
const std::vector< std::string > &  data_names,
const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &  nonscalar_data_ranges,
const DXFlags flags,
std::ostream &  out 
)

Write the given list of patches to the output stream in OpenDX format.

Definition at line 3266 of file data_out_base.cc.

◆ write_eps() [1/2]

template<int spacedim>
void DataOutBase::write_eps ( const std::vector< Patch< 2, spacedim > > &  patches,
const std::vector< std::string > &  data_names,
const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &  nonscalar_data_ranges,
const EpsFlags flags,
std::ostream &  out 
)

Write the given list of patches to the output stream in eps format.

Output in this format circumvents the use of auxiliary graphic programs converting some output format into a graphics format. This has the advantage that output is easy and fast, and the disadvantage that you have to give a whole bunch of parameters which determine the direction of sight, the mode of colorization, the scaling of the height axis, etc. (Of course, all these parameters have reasonable default values, which you may want to change.)

This function only supports output for two-dimensional domains (i.e., with dim==2), with values in the vertical direction taken from a data vector.

Basically, output consists of the mesh and the cells in between them. You can draw either of these, or both, or none if you are really interested in an empty picture. If written, the mesh uses black lines. The cells in between the mesh are either not printed (this will result in a loss of hidden line removal, i.e. you can "see through" the cells to lines behind), printed in white (which does nothing apart from the hidden line removal), or colorized using one of the data vectors (which need not be the same as the one used for computing the height information) and a customizable color function. The default color functions chooses the color between black, blue, green, red and white, with growing values of the data field chosen for colorization. At present, cells are displayed with one color per cell only, which is taken from the value of the data field at the center of the cell; bilinear interpolation of the color on a cell is not used.

By default, the viewpoint is chosen like the default viewpoint in GNUPLOT, i.e. with an angle of 60 degrees with respect to the positive z-axis and rotated 30 degrees in positive sense (as seen from above) away from the negative y-axis. Of course you can change these settings.

EPS output is written without a border around the picture, i.e. the bounding box is close to the output on all four sides. Coordinates are written using at most five digits, to keep picture size at a reasonable size.

All parameters along with their default values are listed in the documentation of the EpsFlags member class of this class. See there for more and detailed information.

Definition at line 4328 of file data_out_base.cc.

◆ write_eps() [2/2]

template<int dim, int spacedim>
void DataOutBase::write_eps ( const std::vector< Patch< dim, spacedim > > &  patches,
const std::vector< std::string > &  data_names,
const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &  nonscalar_data_ranges,
const EpsFlags flags,
std::ostream &  out 
)

This is the same function as above except for domains that are not two-dimensional. This function is not implemented (and will throw an error if called) but is declared to allow for dimension-independent programs.

Definition at line 4310 of file data_out_base.cc.

◆ write_gmv()

template<int dim, int spacedim>
void DataOutBase::write_gmv ( const std::vector< Patch< dim, spacedim > > &  patches,
const std::vector< std::string > &  data_names,
const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &  nonscalar_data_ranges,
const GmvFlags flags,
std::ostream &  out 
)

Write the given list of patches to the output stream in GMV format.

Data is written in the following format: nodes are considered the points of the patches. In spatial dimensions less than three, zeroes are inserted for the missing coordinates. The data vectors are written as node or cell data, where for the first the data space is interpolated to (bi-,tri-)linear elements.

Definition at line 4646 of file data_out_base.cc.

◆ write_gnuplot()

template<int dim, int spacedim>
void DataOutBase::write_gnuplot ( const std::vector< Patch< dim, spacedim > > &  patches,
const std::vector< std::string > &  data_names,
const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &  nonscalar_data_ranges,
const GnuplotFlags flags,
std::ostream &  out 
)

Write the given list of patches to the output stream in gnuplot format. Visualization of two-dimensional data can then be achieved by starting gnuplot and entering the commands

* set data style lines
* splot "filename" using 1:2:n
* 

This example assumes that the number of the data vector displayed is n-2.

The GNUPLOT format is not able to handle data on unstructured grids directly. Directly would mean that you only give the vertices and the solution values thereon and the program constructs its own grid to represent the data. This is only possible for a structured tensor product grid in two dimensions. However, it is possible to give several such patches within one file, which is exactly what the respective function of this class does: writing each cell's data as a patch of data, at least if the patches as passed from derived classes represent cells. Note that the functions on patches need not be continuous at interfaces between patches, so this method also works for discontinuous elements. Note also, that GNUPLOT can do hidden line removal for patched data.

While this discussion applies to two spatial dimensions, it is more complicated in 3d. The reason is that we could still use patches, but it is difficult when trying to visualize them, since if we use a cut through the data (by, for example, using x- and z-coordinates, a fixed y-value and plot function values in z-direction, then the patched data is not a patch in the sense GNUPLOT wants it any more. Therefore, we use another approach, namely writing the data on the 3d grid as a sequence of lines, i.e. two points each associated with one or more data sets. There are therefore 12 lines for each subcells of a patch.

Given the lines as described above, a cut through this data in Gnuplot can then be achieved like this:

*   set data style lines
*   splot [:][:][0:] "T" using 1:2:(\$3==.5 ? \$4 : -1)
* 

This command plots data in \(x\)- and \(y\)-direction unbounded, but in \(z\)-direction only those data points which are above the \(x\)- \(y\)-plane (we assume here a positive solution, if it has negative values, you might want to decrease the lower bound). Furthermore, it only takes the data points with z-values (&3) equal to 0.5, i.e. a cut through the domain at z=0.5. For the data points on this plane, the data values of the first data set (&4) are raised in z-direction above the x-y-plane; all other points are denoted the value -1 instead of the value of the data vector and are not plotted due to the lower bound in z plotting direction, given in the third pair of brackets.

More complex cuts are possible, including nonlinear ones. Note however, that only those points which are actually on the cut-surface are plotted.

Definition at line 3556 of file data_out_base.cc.

◆ write_povray()

template<int dim, int spacedim>
void DataOutBase::write_povray ( const std::vector< Patch< dim, spacedim > > &  patches,
const std::vector< std::string > &  data_names,
const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &  nonscalar_data_ranges,
const PovrayFlags flags,
std::ostream &  out 
)

Write the given list of patches to the output stream for the Povray raytracer.

Output in this format creates a povray source file, include standard camera and light source definition for rendering with povray 3.1 At present, this format only supports output for two-dimensional data, with values in the third direction taken from a data vector.

The output uses two different povray-objects:

  • BICUBIC_PATCH A bicubic_patch is a 3-dimensional Bezier patch. It consists of 16 Points describing the surface. The 4 corner points are touched by the object, while the other 12 points pull and stretch the patch into shape. One bicubic_patch is generated on each patch. Therefore the number of subdivisions has to be 3 to provide the patch with 16 points. A bicubic patch is not exact but generates very smooth images.

  • MESH The mesh object is used to store large number of triangles. Every square of the patch data is split into one upper-left and one lower-right triangle. If the number of subdivisions is three, 32 triangle are generated for every patch.

    Using the smooth flag povray interpolates the normals on the triangles, imitating a curved surface

All objects get one texture definition called Tex. This texture has to be declared somewhere before the object data. This may be in an external data file or at the beginning of the output file. Setting the external_data flag to false, an standard camera, light and texture (scaled to fit the scene) is added to the output file. Set to true an include file "data.inc" is included. This file is not generated by deal and has to include camera, light and the texture definition Tex.

You need povray (>=3.0) to render the scene. The minimum options for povray are:

*   povray +I<inputfile> +W<horiz. size> +H<ver. size> +L<include path>
* 

If the external file "data.inc" is used, the path to this file has to be included in the povray options.

Definition at line 4292 of file data_out_base.cc.

◆ write_tecplot()

template<int dim, int spacedim>
void DataOutBase::write_tecplot ( const std::vector< Patch< dim, spacedim > > &  patches,
const std::vector< std::string > &  data_names,
const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &  nonscalar_data_ranges,
const TecplotFlags flags,
std::ostream &  out 
)

Write the given list of patches to the output stream in Tecplot ASCII format (FEBLOCK).

For more information consult the Tecplot Users and Reference manuals.

Definition at line 4779 of file data_out_base.cc.

◆ write_ucd()

template<int dim, int spacedim>
void DataOutBase::write_ucd ( const std::vector< Patch< dim, spacedim > > &  patches,
const std::vector< std::string > &  data_names,
const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &  nonscalar_data_ranges,
const UcdFlags flags,
std::ostream &  out 
)

Write the given list of patches to the output stream in UCD format described in the AVS developer's guide (now AVS). Due to limitations in the present format, only node based data can be output, which in one reason why we invented the patch concept. In order to write higher order elements, you may split them up into several subdivisions of each cell. These subcells will then, however, also appear as different cells by programs which understand the UCD format.

No use is made of the possibility to give model data since these are not supported by all UCD aware programs. You may give cell data in derived classes by setting all values of a given data set on a patch to the same value.

Definition at line 3175 of file data_out_base.cc.

◆ write_vtk()

template<int dim, int spacedim>
void DataOutBase::write_vtk ( const std::vector< Patch< dim, spacedim > > &  patches,
const std::vector< std::string > &  data_names,
const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &  nonscalar_data_ranges,
const VtkFlags flags,
std::ostream &  out 
)

Write the given list of patches to the output stream in VTK format. The data is written in the traditional VTK format as opposed to the XML-based format that write_vtu() produces.

The nonscalar_data_ranges argument denotes ranges of components in the output that are considered a vector, rather than simply a collection of scalar fields. The VTK output format has special provisions that allow these components to be output by a single name rather than having to group several scalar fields into a vector later on in the visualization program.

Note
VTK is a legacy format and has largely been supplanted by the VTU format (an XML-structured version of VTK). In particular, VTU allows for the compression of data and consequently leads to much smaller file sizes that equivalent VTK files for large files. Since all visualization programs that support VTK also support VTU, you should consider using the latter file format instead, by using the write_vtu() function.

Definition at line 4932 of file data_out_base.cc.

◆ write_vtu()

template<int dim, int spacedim>
void DataOutBase::write_vtu ( const std::vector< Patch< dim, spacedim > > &  patches,
const std::vector< std::string > &  data_names,
const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &  nonscalar_data_ranges,
const VtkFlags flags,
std::ostream &  out 
)

Write the given list of patches to the output stream in VTU format. The data is written in the XML-based VTK format as opposed to the traditional format that write_vtk() produces.

The nonscalar_data_ranges argument denotes ranges of components in the output that are considered a vector, rather than simply a collection of scalar fields. The VTK output format has special provisions that allow these components to be output by a single name rather than having to group several scalar fields into a vector later on in the visualization program.

Some visualization programs, such as ParaView, can read several separate VTU files to parallelize visualization. In that case, you need a .pvtu file that describes which VTU files form a group. The DataOutInterface::write_pvtu_record() function can generate such a centralized record. Likewise, DataOutInterface::write_visit_record() does the same for VisIt (although VisIt can also read pvtu records since version 2.5.1). Finally, for time dependent problems, you may also want to look at DataOutInterface::write_pvd_record()

The use of this function is explained in step-40.

Definition at line 5228 of file data_out_base.cc.

◆ write_vtu_header()

void DataOutBase::write_vtu_header ( std::ostream &  out,
const VtkFlags flags 
)

This writes the header for the xml based vtu file format. This routine is used internally together with DataOutInterface::write_vtu_footer() and DataOutInterface::write_vtu_main() by DataOutBase::write_vtu().

Definition at line 5180 of file data_out_base.cc.

◆ write_vtu_footer()

void DataOutBase::write_vtu_footer ( std::ostream &  out)

This function writes the footer for the xml based vtu file format. This routine is used internally together with DataOutInterface::write_vtu_header() and DataOutInterface::write_vtu_main() by DataOutBase::write_vtu().

Definition at line 5217 of file data_out_base.cc.

◆ write_vtu_main()

template<int dim, int spacedim>
void DataOutBase::write_vtu_main ( const std::vector< Patch< dim, spacedim > > &  patches,
const std::vector< std::string > &  data_names,
const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &  nonscalar_data_ranges,
const VtkFlags flags,
std::ostream &  out 
)

This function writes the main part for the xml based vtu file format. This routine is used internally together with DataOutInterface::write_vtu_header() and DataOutInterface::write_vtu_footer() by DataOutBase::write_vtu().

Definition at line 5250 of file data_out_base.cc.

◆ write_pvtu_record()

void DataOutBase::write_pvtu_record ( std::ostream &  out,
const std::vector< std::string > &  piece_names,
const std::vector< std::string > &  data_names,
const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &  nonscalar_data_ranges,
const VtkFlags flags 
)

Some visualization programs, such as ParaView, can read several separate VTU files that all form part of the same simulation, in order to parallelize visualization. In that case, you need a .pvtu file that describes which VTU files (written, for example, through the DataOutInterface::write_vtu() function) form a group. The current function can generate such a centralized record.

This function is typically not called by itself from user space, but you may want to call it through DataOutInterface::write_pvtu_record() since the DataOutInterface class has access to information that you would have to provide to the current function by hand.

In any case, whether this function is called directly or via DataOutInterface::write_pvtu_record(), the central record file so written contains a list of (scalar or vector) fields that describes which fields can actually be found in the individual files that comprise the set of parallel VTU files along with the names of these files. This function gets the names and types of fields through the third and fourth argument; you can determine these by hand, but in practice, this function is most easily called by calling DataOutInterfaces::write_pvtu_record(), which determines the last two arguments by calling DataOutInterface::get_dataset_names() and DataOutInterface::get_nonscalar_data_ranges() functions. The second argument to this function specifies the names of the files that form the parallel set.

Note
Use DataOutBase::write_vtu() and DataOutInterface::write_vtu() for writing each piece. Also note that only one parallel process needs to call the current function, listing the names of the files written by all parallel processes.
In order to tell Paraview to group together multiple pvtu files that each describe one time step of a time dependent simulation, see the DataOutBase::write_pvd_record() function.
Older versions of VisIt (before 2.5.1), can not read pvtu records. However, it can read visit records as written by the write_visit_record() function.

Definition at line 6103 of file data_out_base.cc.

◆ write_pvd_record()

void DataOutBase::write_pvd_record ( std::ostream &  out,
const std::vector< std::pair< double, std::string > > &  times_and_names 
)

In ParaView it is possible to visualize time-dependent data tagged with the current integration time of a time dependent simulation. To use this feature you need a .pvd file that describes which VTU or PVTU file belongs to which timestep. This function writes a file that provides this mapping, i.e., it takes a list of pairs each of which indicates a particular time instant and the corresponding file that contains the graphical data for this time instant.

A typical use case, in program that computes a time dependent solution, would be the following (time and time_step are member variables of the class with types double and unsigned int, respectively; the variable times_and_names is of type std::vector<std::pair<double,std::string> >):

template <int dim>
void MyEquation<dim>::output_results () const
{
DataOut<dim> data_out;
data_out.attach_dof_handler(dof_handler);
data_out.add_data_vector(solution, "U");
data_out.build_patches();
const std::string filename = "solution-" +
Utilities::int_to_string (timestep_n, 3) +
".vtu";
std::ofstream output(filename);
data_out.write_vtu(output);
times_and_names.emplace_back (time, filename);
std::ofstream pvd_output ("solution.pvd");
DataOutBase::write_pvd_record (pvd_output, times_and_names);
}
void write_vtu(std::ostream &out) const
void attach_dof_handler(const DoFHandler< dim, spacedim > &)
void add_data_vector(const VectorType &data, const std::vector< std::string > &names, const DataVectorType type=type_automatic, const std::vector< DataComponentInterpretation::DataComponentInterpretation > &data_component_interpretation={})
virtual void build_patches(const unsigned int n_subdivisions=0)
Definition data_out.cc:1062
void write_pvd_record(std::ostream &out, const std::vector< std::pair< double, std::string > > &times_and_names)
std::string int_to_string(const unsigned int value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition utilities.cc:470
Note
See DataOutInterface::write_vtu, DataOutInterface::write_pvtu_record, and DataOutInterface::write_vtu_in_parallel for writing solutions at each timestep.
The second element of each pair, i.e., the file in which the graphical data for each time is stored, may itself be again a file that references other files. For example, it could be the name for a .pvtu file that references multiple parts of a parallel computation.

Definition at line 6256 of file data_out_base.cc.

◆ write_visit_record() [1/3]

void DataOutBase::write_visit_record ( std::ostream &  out,
const std::vector< std::string > &  piece_names 
)

This function is the exact equivalent of the write_pvtu_record() function but for older versions of the VisIt visualization program and for one visualization graph (or one time step only). See there for the purpose of this function.

This function is documented in the "Creating a master file for parallel" section (section 5.7) of the "Getting data into VisIt" report that can be found here: https://wci.llnl.gov/codes/visit/2.0.0/GettingDataIntoVisIt2.0.0.pdf

Definition at line 6293 of file data_out_base.cc.

◆ write_visit_record() [2/3]

void DataOutBase::write_visit_record ( std::ostream &  out,
const std::vector< std::vector< std::string > > &  piece_names 
)

This function is equivalent to the write_visit_record() above but for multiple time steps. Here is an example of how the function would be used:

const unsigned int number_of_time_steps = 3;
std::vector<std::vector<std::string > > piece_names(number_of_time_steps);
piece_names[0].emplace_back("subdomain_01.time_step_0.vtk");
piece_names[0].emplace_back("subdomain_02.time_step_0.vtk");
piece_names[1].emplace_back("subdomain_01.time_step_1.vtk");
piece_names[1].emplace_back("subdomain_02.time_step_1.vtk");
piece_names[2].emplace_back("subdomain_01.time_step_2.vtk");
piece_names[2].emplace_back("subdomain_02.time_step_2.vtk");
std::ofstream visit_output ("solution.visit");
DataOutBase::write_visit_record(visit_output, piece_names);
void write_visit_record(std::ostream &out, const std::vector< std::string > &piece_names)

This function is documented in the "Creating a master file for parallel" section (section 5.7) of the "Getting data into VisIt" report that can be found here: https://wci.llnl.gov/codes/visit/2.0.0/GettingDataIntoVisIt2.0.0.pdf

Definition at line 6306 of file data_out_base.cc.

◆ write_visit_record() [3/3]

void DataOutBase::write_visit_record ( std::ostream &  out,
const std::vector< std::pair< double, std::vector< std::string > > > &  times_and_piece_names 
)

This function is equivalent to the write_visit_record() above but for multiple time steps and with additional information about the time for each timestep. Here is an example of how the function would be used:

const unsigned int number_of_time_steps = 3;
std::vector<std::pair<double,std::vector<std::string > > >
times_and_piece_names(number_of_time_steps);
times_and_piece_names[0].first = 0.0;
times_and_piece_names[0].second.emplace_back("subdomain_01.time_step_0.vtk");
times_and_piece_names[0].second.emplace_back("subdomain_02.time_step_0.vtk");
times_and_piece_names[1].first = 0.5;
times_and_piece_names[1].second.emplace_back("subdomain_01.time_step_1.vtk");
times_and_piece_names[1].second.emplace_back("subdomain_02.time_step_1.vtk");
times_and_piece_names[2].first = 1.0;
times_and_piece_names[2].second.emplace_back("subdomain_01.time_step_2.vtk");
times_and_piece_names[2].second.emplace_back("subdomain_02.time_step_2.vtk");
std::ofstream visit_output ("solution.visit");
DataOutBase::write_visit_record(visit_output, times_and_piece_names);

This function is documented in the "Creating a master file for parallel" section (section 5.7) of the "Getting data into VisIt" report that can be found here: https://wci.llnl.gov/codes/visit/2.0.0/GettingDataIntoVisIt2.0.0.pdf

Definition at line 6334 of file data_out_base.cc.

◆ write_svg() [1/2]

template<int spacedim>
void DataOutBase::write_svg ( const std::vector< Patch< 2, spacedim > > &  patches,
const std::vector< std::string > &  data_names,
const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &  nonscalar_data_ranges,
const SvgFlags flags,
std::ostream &  out 
)

Write the given list of patches to the output stream in SVG format.

SVG (Scalable Vector Graphics) is an XML-based vector image format developed and maintained by the World Wide Web Consortium (W3C). This function conforms to the latest specification SVG 1.1, released on August 16, 2011. Controlling the graphic output is possible by setting or clearing the respective flags (see the SvgFlags struct). At present, this format only supports output for two-dimensional data, with values in the third direction taken from a data vector.

For the output, each patch is subdivided into four triangles which are then written as polygons and filled with a linear color gradient. The arising coloring of the patches visualizes the data values at the vertices taken from the specified data vector. A colorbar can be drawn to encode the coloring.

Note
This function is so far only implemented for two dimensions with an additional dimension reserved for data information.

Definition at line 6386 of file data_out_base.cc.

◆ write_deal_II_intermediate()

template<int dim, int spacedim>
void DataOutBase::write_deal_II_intermediate ( const std::vector< Patch< dim, spacedim > > &  patches,
const std::vector< std::string > &  data_names,
const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &  nonscalar_data_ranges,
const Deal_II_IntermediateFlags flags,
std::ostream &  out 
)

Write the given list of patches to the output stream in deal.II intermediate format. This is not a format understood by any other graphics program, but is rather a direct dump of the intermediate internal format used by deal.II. This internal format is generated by the various classes that can generate output using the DataOutBase class, for example from a finite element solution, and is then converted in the present class to the final graphics format.

Note that the intermediate format is what its name suggests: a direct representation of internal data. It isn't standardized and will change whenever we change our internal representation. You can only expect to process files written in this format using the same version of deal.II that was used for writing.

The reason why we offer to write out this intermediate format is that it can be read back into a deal.II program using the DataOutReader class, which is helpful in at least two contexts: First, this can be used to later generate graphical output in any other graphics format presently understood; this way, it is not necessary to know at run-time which output format is requested, or if multiple output files in different formats are needed. Secondly, in contrast to almost all other graphics formats, it is possible to merge several files that contain intermediate format data, and generate a single output file from it, which may be again in intermediate format or any of the final formats. This latter option is most helpful for parallel programs: as demonstrated in the step-17 example program, it is possible to let only one processor generate the graphical output for the entire parallel program, but this can become vastly inefficient if many processors are involved, because the load is no longer balanced. The way out is to let each processor generate intermediate graphical output for its chunk of the domain, and the later merge the different files into one, which is an operation that is much cheaper than the generation of the intermediate data.

Intermediate format deal.II data is usually stored in files with the ending .d2.

Definition at line 7352 of file data_out_base.cc.

◆ write_deal_II_intermediate_in_parallel()

template<int dim, int spacedim>
void DataOutBase::write_deal_II_intermediate_in_parallel ( const std::vector< Patch< dim, spacedim > > &  patches,
const std::vector< std::string > &  data_names,
const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &  nonscalar_data_ranges,
const Deal_II_IntermediateFlags flags,
const std::string &  filename,
const MPI_Comm  comm,
const CompressionLevel  compression 
)

Like write_deal_II_intermediate() but write all patches from all ranks using MPI I/O into a single file with name name. Compression using zlib is optional and controlled by the compression argument.

The files typically have the extension .pd2.

Definition at line 7400 of file data_out_base.cc.

◆ write_hdf5_parallel() [1/2]

template<int dim, int spacedim>
void DataOutBase::write_hdf5_parallel ( const std::vector< Patch< dim, spacedim > > &  patches,
const DataOutFilter data_filter,
const DataOutBase::Hdf5Flags flags,
const std::string &  filename,
const MPI_Comm  comm 
)

Write the data in data_filter to a single HDF5 file containing both the mesh and solution values.

Definition at line 8689 of file data_out_base.cc.

◆ write_hdf5_parallel() [2/2]

template<int dim, int spacedim>
void DataOutBase::write_hdf5_parallel ( const std::vector< Patch< dim, spacedim > > &  patches,
const DataOutFilter data_filter,
const DataOutBase::Hdf5Flags flags,
const bool  write_mesh_file,
const std::string &  mesh_filename,
const std::string &  solution_filename,
const MPI_Comm  comm 
)

Write the data in data_filter to HDF5 file(s). If write_mesh_file is false, the mesh data will not be written and the solution file will contain only the solution values. If write_mesh_file is true and the filenames are the same, the resulting file will contain both mesh data and solution values.

Definition at line 8704 of file data_out_base.cc.

◆ write_filtered_data()

template<int dim, int spacedim>
void DataOutBase::write_filtered_data ( const std::vector< Patch< dim, spacedim > > &  patches,
const std::vector< std::string > &  data_names,
const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &  nonscalar_data_ranges,
DataOutBase::DataOutFilter filtered_data 
)

DataOutFilter is an intermediate data format that reduces the amount of data that will be written to files. The object filled by this function can then later be used again to write data in a concrete file format; see, for example, DataOutBase::write_hdf5_parallel().

Definition at line 8533 of file data_out_base.cc.

◆ determine_intermediate_format_dimensions()

std::pair< unsigned int, unsigned int > DataOutBase::determine_intermediate_format_dimensions ( std::istream &  input)

Given an input stream that contains data written by write_deal_II_intermediate(), determine the dim and spacedim template parameters with which that function was called, and return them as a pair of values.

Note that this function eats a number of elements at the present position of the stream, and therefore alters it. In order to read from it using, for example, the DataOutReader class, you may wish to either reset the stream to its previous position, or close and reopen it.

Definition at line 7564 of file data_out_base.cc.

◆ parse_output_format()

OutputFormat DataOutBase::parse_output_format ( const std::string &  format_name)

Return the OutputFormat value corresponding to the given string. If the string does not match any known format, an exception is thrown.

The main purpose of this function is to allow a program to use any implemented output format without the need to extend the program's parser each time a new format is implemented.

To get a list of presently available format names, e.g. to give it to the ParameterHandler class, use the function get_output_format_names().

Definition at line 2439 of file data_out_base.cc.

◆ get_output_format_names()

std::string DataOutBase::get_output_format_names ( )

Return a list of implemented output formats. The different names are separated by vertical bar signs (‘|’) as used by the ParameterHandler classes.

Definition at line 2488 of file data_out_base.cc.

◆ default_suffix()

std::string DataOutBase::default_suffix ( const OutputFormat  output_format)

Provide a function which tells us which suffix a file with a given output format usually has. At present the following formats are defined:

  • dx: .dx
  • ucd: .inp
  • gnuplot: .gnuplot
  • povray: .pov
  • eps: .eps
  • gmv: .gmv
  • tecplot: .dat
  • vtk: .vtk
  • vtu: .vtu
  • svg: .svg
  • deal_II_intermediate: .d2.
Deprecated:
Using Tecplot binary output is deprecated.

Definition at line 2496 of file data_out_base.cc.

◆ operator<<()

template<int dim, int spacedim>
std::ostream & DataOutBase::operator<< ( std::ostream &  out,
const Patch< dim, spacedim > &  patch 
)
private

Output operator for an object of type DataOutBase::Patch. This operator dumps the intermediate graphics format represented by the patch data structure. It may later be converted into regular formats for a number of graphics programs.

Definition at line 9684 of file data_out_base.cc.

◆ operator>>()

template<int dim, int spacedim>
std::istream & DataOutBase::operator>> ( std::istream &  in,
Patch< dim, spacedim > &  patch 
)
private

Input operator for an object of type DataOutBase::Patch. This operator reads the intermediate graphics format represented by the patch data structure, using the format in which it was written using the operator<<.

Definition at line 9721 of file data_out_base.cc.

◆ get_node_positions()

template<int dim, int spacedim>
std::vector< Point< spacedim > > DataOutBase::get_node_positions ( const std::vector< Patch< dim, spacedim > > &  patches)

Obtain the positions of all nodes referenced in the patches given as argument.

Definition at line 2542 of file data_out_base.cc.

◆ write_nodes()

template<int dim, int spacedim, typename StreamType >
void DataOutBase::write_nodes ( const std::vector< Patch< dim, spacedim > > &  patches,
StreamType &  out 
)

Definition at line 2603 of file data_out_base.cc.

◆ write_cells()

template<int dim, int spacedim, typename StreamType >
void DataOutBase::write_cells ( const std::vector< Patch< dim, spacedim > > &  patches,
StreamType &  out 
)

Definition at line 2620 of file data_out_base.cc.

◆ write_high_order_cells()

template<int dim, int spacedim, typename StreamType >
void DataOutBase::write_high_order_cells ( const std::vector< Patch< dim, spacedim > > &  patches,
StreamType &  out,
const bool  legacy_format 
)

Definition at line 2719 of file data_out_base.cc.

◆ write_data()

template<int dim, int spacedim, typename StreamType >
void DataOutBase::write_data ( const std::vector< Patch< dim, spacedim > > &  patches,
unsigned int  n_data_sets,
const bool  double_precision,
StreamType &  out 
)

Definition at line 2830 of file data_out_base.cc.

◆ write_svg() [2/2]

template<int dim, int spacedim>
void DataOutBase::write_svg ( const std::vector< Patch< dim, spacedim > > &  ,
const std::vector< std::string > &  ,
const std::vector< std::tuple< unsigned int, unsigned int, std::string, DataComponentInterpretation::DataComponentInterpretation > > &  ,
const SvgFlags ,
std::ostream &   
)

Definition at line 6370 of file data_out_base.cc.