* (bi-,tri-)linear elements.
*
*
+ * @sect3{VTK format}
+ *
+ * This is the file format used by the Visualization Toolkit (VTK, see
+ * @url{http://www.kitware.com/vtk.html}), as
+ * described in their manual, section 14.3. It is similar to the GMV
+ * format, see there for more information.
+ *
+ *
* @sect2{Output parameters}
*
* All functions take a parameter which is a structure of type @p{XFlags}, where
* present no flags are implemented.
*/
struct GmvFlags
+ {
+ private:
+ /**
+ * Dummy entry to suppress compiler
+ * warnings when copying an empty
+ * structure. Remove this member
+ * when adding the first flag to
+ * this structure (and remove the
+ * @p{private} as well).
+ */
+ int dummy;
+
+ public:
+ /**
+ * Declare all flags with name
+ * and type as offered by this
+ * class, for use in input files.
+ */
+ static void declare_parameters (ParameterHandler &prm);
+
+ /**
+ * Read the parameters declared in
+ * @p{declare_parameters} and set the
+ * flags for this output format
+ * accordingly.
+ *
+ * The flags thus obtained overwrite
+ * all previous contents of this object.
+ */
+ void parse_parameters (ParameterHandler &prm);
+
+ /**
+ * Determine an estimate for
+ * the memory consumption (in
+ * bytes) of this
+ * object. Since sometimes
+ * the size of objects can
+ * not be determined exactly
+ * (for example: what is the
+ * memory consumption of an
+ * STL @p{std::map} type with a
+ * certain number of
+ * elements?), this is only
+ * an estimate. however often
+ * quite close to the true
+ * value.
+ */
+ unsigned int memory_consumption () const;
+ };
+
+ /**
+ * Flags describing the details
+ * of output in VTK format. At
+ * present no flags are
+ * implemented.
+ */
+ struct VtkFlags
{
private:
/**
static void write_gmv (const typename std::vector<Patch<dim,spacedim> > &patches,
const std::vector<std::string> &data_names,
const GmvFlags &flags,
+ std::ostream &out);
+
+ /**
+ * Write the given list of
+ * patches to the output stream
+ * in vtk format. See the general
+ * documentation for more
+ * information on the parameters.
+ */
+ template <int dim, int spacedim>
+ static void write_vtk (const typename std::vector<Patch<dim,spacedim> > &patches,
+ const std::vector<std::string> &data_names,
+ const VtkFlags &flags,
std::ostream &out);
/**
* coordinates, one height value)
* to the direction of sight.
*/
- class EpsCell2d {
+ class EpsCell2d
+ {
public:
/**
* from the main thread, and is
* thus moved into this separate
* function.
+ *
+ * Note that because of the close
+ * similarity of the two formats,
+ * this function is also used by
+ * the Vtk output function.
*/
template <int dim, int spacedim>
static void
{
public:
/**
- * Provide a data type specifying the
- * presently supported output formats.
+ * Provide a data type specifying
+ * the presently supported output
+ * formats.
*/
- enum OutputFormat { default_format, ucd, gnuplot, povray, eps, gmv };
+ enum OutputFormat { default_format, ucd, gnuplot, povray, eps, gmv, vtk };
/**
- * Obtain data through the @p{get_patches}
- * function and write it to @p{out} in
- * UCD format.
+ * Obtain data through the
+ * @p{get_patches} function and
+ * write it to @p{out} in UCD
+ * format.
*/
void write_ucd (std::ostream &out) const;
/**
- * Obtain data through the @p{get_patches}
- * function and write it to @p{out} in
- * GNUPLOT format.
+ * Obtain data through the
+ * @p{get_patches} function and
+ * write it to @p{out} in GNUPLOT
+ * format.
*/
void write_gnuplot (std::ostream &out) const;
/**
- * Obtain data through the @p{get_patches}
- * function and write it to @p{out} in
- * POVRAY format.
+ * Obtain data through the
+ * @p{get_patches} function and
+ * write it to @p{out} in POVRAY
+ * format.
*/
void write_povray (std::ostream &out) const;
/**
- * Obtain data through the @p{get_patches}
- * function and write it to @p{out} in
- * EPS format.
+ * Obtain data through the
+ * @p{get_patches} function and
+ * write it to @p{out} in EPS
+ * format.
*/
void write_eps (std::ostream &out) const;
/**
- * Obtain data through the @p{get_patches}
- * function and write it to @p{out} in
- * GMV format.
+ * Obtain data through the
+ * @p{get_patches} function and
+ * write it to @p{out} in GMV
+ * format.
*/
void write_gmv (std::ostream &out) const;
+ /**
+ * Obtain data through the
+ * @p{get_patches} function and
+ * write it to @p{out} in Vtk
+ * format.
+ */
+ void write_vtk (std::ostream &out) const;
+
/**
- * Write data and grid to @p{out} according
- * to the given data format. This function
- * simply calls the appropriate
- * @p{write_*} function. If no output format is
- * requested, the @p{default_format} is written.
+ * Write data and grid to @p{out}
+ * according to the given data
+ * format. This function simply
+ * calls the appropriate
+ * @p{write_*} function. If no
+ * output format is requested,
+ * the @p{default_format} is
+ * written.
*
- * An error occurs if no format is provided and
- * the default format is @p{default_format}.
+ * An error occurs if no format
+ * is provided and the default
+ * format is @p{default_format}.
*/
void write (std::ostream &out,
const OutputFormat output_format = default_format) const;
/**
- * Set the default format. The value set here
- * is used anytime, output for format @p{default_format} is
+ * Set the default format. The
+ * value set here is used
+ * anytime, output for format
+ * @p{default_format} is
* requested.
*/
- void set_default_format(const OutputFormat default_format);
+ void set_default_format (const OutputFormat default_format);
/**
- * Set the flags to be used for output
- * in UCD format.
+ * Set the flags to be used for
+ * output in UCD format.
*/
void set_flags (const UcdFlags &ucd_flags);
/**
- * Set the flags to be used for output
- * in GNUPLOT format.
+ * Set the flags to be used for
+ * output in GNUPLOT format.
*/
void set_flags (const GnuplotFlags &gnuplot_flags);
/**
- * Set the flags to be used for output
- * in POVRAY format.
+ * Set the flags to be used for
+ * output in POVRAY format.
*/
void set_flags (const PovrayFlags &povray_flags);
/**
- * Set the flags to be used for output
- * in EPS output.
+ * Set the flags to be used for
+ * output in EPS output.
*/
void set_flags (const EpsFlags &eps_flags);
/**
- * Set the flags to be used for output
- * in GMV format.
+ * Set the flags to be used for
+ * output in GMV format.
*/
void set_flags (const GmvFlags &gmv_flags);
+ /**
+ * Set the flags to be used for
+ * output in GMV format.
+ */
+ void set_flags (const VtkFlags &vtk_flags);
+
/**
* Provide a function which tells us which
* @item @p{gnuplot}: @p{.gnuplot}
* @item @p{povray}: @p{.pov}
* @item @p{eps}: @p{.eps}
- * @item @p{gmv}: @p{.gmv}.
+ * @item @p{gmv}: @p{.gmv}
+ * @item @p{vtk}: @p{.vtk}.
* @end{itemize}
*
* If this function is called
- * with no argument or @p{default_format}, the
- * suffix for the
- * @p{default_format} is returned.
+ * with no argument or
+ * @p{default_format}, the suffix
+ * for the @p{default_format} is
+ * returned.
*/
std::string default_suffix (const OutputFormat output_format = default_format) const;
/**
- * Return the @p{OutputFormat} value
- * corresponding to the given string. If
- * the string does not match any known
- * format, an exception is thrown.
+ * Return the @p{OutputFormat}
+ * value corresponding to the
+ * given string. If the string
+ * does not match any known
+ * format, an exception is
+ * thrown.
*
- * Since this function does not need data
- * from this object, it is static and can
- * thus be called without creating an
- * object of this class. Its main purpose
- * 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.
+ * Since this function does not
+ * need data from this object, it
+ * is static and can thus be
+ * called without creating an
+ * object of this class. Its main
+ * purpose 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
- * @p{ParameterHandler} class, use the
- * function @p{get_output_format_names ()}.
+ * To get a list of presently
+ * available format names,
+ * e.g. to give it to the
+ * @p{ParameterHandler} class,
+ * use the function
+ * @p{get_output_format_names
+ * ()}.
*/
static OutputFormat parse_output_format (const std::string &format_name);
/**
- * Return a list of implemented output
- * formats. The different names are
- * separated by vertical bar signs (@p{`|'})
- * as used by the @p{ParameterHandler}
- * classes.
+ * Return a list of implemented
+ * output formats. The different
+ * names are separated by
+ * vertical bar signs (@p{`|'})
+ * as used by the
+ * @p{ParameterHandler} classes.
*/
static std::string get_output_format_names ();
/**
- * Declare parameters for all output
- * formats by declaring subsections
- * within the parameter file for each
- * output format and call the
- * respective @p{declare_parameters}
- * functions of the flag classes for
- * each output format.
+ * Declare parameters for all
+ * output formats by declaring
+ * subsections within the
+ * parameter file for each output
+ * format and call the respective
+ * @p{declare_parameters}
+ * functions of the flag classes
+ * for each output format.
*
- * Some of the declared subsections
- * may not contain entries, if the
- * respective format does not
- * export any flags.
+ * Some of the declared
+ * subsections may not contain
+ * entries, if the respective
+ * format does not export any
+ * flags.
*
- * Note that the top-level parameters
- * denoting the number of subdivisions
- * per patch and the output format
- * are not declared, since they are
- * only passed to virtual functions
- * and are not stored inside objects
- * of this type. You have to
- * declare them yourself.
+ * Note that the top-level
+ * parameters denoting the number
+ * of subdivisions per patch and
+ * the output format are not
+ * declared, since they are only
+ * passed to virtual functions
+ * and are not stored inside
+ * objects of this type. You have
+ * to declare them yourself.
*/
static void declare_parameters (ParameterHandler &prm);
/**
- * Read the parameters declared in
- * @p{declare_parameters} and set the
- * flags for the output formats
- * accordingly.
+ * Read the parameters declared
+ * in @p{declare_parameters} and
+ * set the flags for the output
+ * formats accordingly.
*
- * The flags thus obtained overwrite
- * all previous contents of the flag
- * objects as default-constructed or
- * set by the @p{set_flags} function.
+ * The flags thus obtained
+ * overwrite all previous
+ * contents of the flag objects
+ * as default-constructed or set
+ * by the @p{set_flags} function.
*/
void parse_parameters (ParameterHandler &prm);
protected:
/**
- * This is the abstract function through
- * which derived classes propagate
- * preprocessed data in the form of
- * @p{Patch} structures (declared in
- * the base class @p{DataOutBase}) to
- * the actual output function. You
- * need to overload this function to
- * allow the output functions to
- * know what they shall print.
+ * This is the abstract function
+ * through which derived classes
+ * propagate preprocessed data in
+ * the form of @p{Patch}
+ * structures (declared in the
+ * base class @p{DataOutBase}) to
+ * the actual output
+ * function. You need to overload
+ * this function to allow the
+ * output functions to know what
+ * they shall print.
*/
virtual const typename std::vector<typename DataOutBase::Patch<dim,spacedim> > &
get_patches () const = 0;
/**
- * Abstract virtual function through
- * which the names of data sets are
- * obtained by the output functions
- * of the base class.
+ * Abstract virtual function
+ * through which the names of
+ * data sets are obtained by the
+ * output functions of the base
+ * class.
*/
virtual std::vector<std::string> get_dataset_names () const = 0;
OutputFormat default_fmt;
/**
- * Flags to be used upon output of UCD
- * data. Can be changed by using the
- * @p{set_flags} function.
+ * Flags to be used upon output
+ * of UCD data. Can be changed by
+ * using the @p{set_flags}
+ * function.
*/
UcdFlags ucd_flags;
/**
- * Flags to be used upon output of GNUPLOT
- * data. Can be changed by using the
+ * Flags to be used upon output
+ * of GNUPLOT data. Can be
+ * changed by using the
* @p{set_flags} function.
*/
GnuplotFlags gnuplot_flags;
/**
- * Flags to be used upon output of POVRAY
- * data. Can be changed by using the
- * @p{set_flags} function.
+ * Flags to be used upon output
+ * of POVRAY data. Can be changed
+ * by using the @p{set_flags}
+ * function.
*/
PovrayFlags povray_flags;
/**
- * Flags to be used upon output of EPS
- * data in one space dimension. Can be
- * changed by using the @p{set_flags}
+ * Flags to be used upon output
+ * of EPS data in one space
+ * dimension. Can be changed by
+ * using the @p{set_flags}
* function.
*/
EpsFlags eps_flags;
/**
- * Flags to be used upon output of gmv
- * data in one space dimension. Can be
- * changed by using the @p{set_flags}
+ * Flags to be used upon output
+ * of gmv data in one space
+ * dimension. Can be changed by
+ * using the @p{set_flags}
* function.
*/
GmvFlags gmv_flags;
+
+ /**
+ * Flags to be used upon output
+ * of vtk data in one space
+ * dimension. Can be changed by
+ * using the @p{set_flags}
+ * function.
+ */
+ VtkFlags vtk_flags;
};
+template <int dim, int spacedim>
+void DataOutBase::write_vtk (const typename std::vector<Patch<dim,spacedim> > &patches,
+ const std::vector<std::string> &data_names,
+ const VtkFlags &/*flags*/,
+ std::ostream &out)
+{
+ AssertThrow (out, ExcIO());
+
+ Assert (patches.size() > 0, ExcNoPatches());
+
+ const unsigned int n_data_sets = data_names.size();
+ // check against # of data sets in
+ // first patch. checks against all
+ // other patches are made in
+ // write_gmv_reorder_data_vectors
+ Assert (n_data_sets == patches[0].data.m(),
+ ExcUnexpectedNumberOfDatasets (patches[0].data.m(), n_data_sets));
+
+
+ ///////////////////////
+ // preamble
+ if (true)
+ {
+ time_t time1= time (0);
+ tm *time = localtime(&time1);
+ out << "# vtk DataFile Version 3.0"
+ << std::endl
+ << "This file was generated by the deal.II library on "
+ << time->tm_year+1900 << "/"
+ << time->tm_mon+1 << "/"
+ << time->tm_mday << " at "
+ << time->tm_hour << ":"
+ << std::setw(2) << time->tm_min << ":"
+ << std::setw(2) << time->tm_sec
+ << std::endl
+ << "ASCII"
+ << std::endl
+ << "DATASET UNSTRUCTURED_GRID"
+ << std::endl;
+ };
+
+
+ // first count the number of cells
+ // and cells for later use
+ unsigned int n_cells = 0,
+ n_nodes = 0;
+ for (typename std::vector<Patch<dim,spacedim> >::const_iterator patch=patches.begin();
+ patch!=patches.end(); ++patch)
+ switch (dim)
+ {
+ case 1:
+ n_cells += patch->n_subdivisions;
+ n_nodes += patch->n_subdivisions+1;
+ break;
+ case 2:
+ n_cells += patch->n_subdivisions *
+ patch->n_subdivisions;
+ n_nodes += (patch->n_subdivisions+1) *
+ (patch->n_subdivisions+1);
+ break;
+ case 3:
+ n_cells += patch->n_subdivisions *
+ patch->n_subdivisions *
+ patch->n_subdivisions;
+ n_nodes += (patch->n_subdivisions+1) *
+ (patch->n_subdivisions+1) *
+ (patch->n_subdivisions+1);
+ break;
+ default:
+ Assert (false, ExcNotImplemented());
+ };
+
+
+ // in gmv format the vertex
+ // coordinates and the data have an
+ // order that is a bit unpleasant
+ // (first all x coordinates, then
+ // all y coordinate, ...; first all
+ // data of variable 1, then
+ // variable 2, etc), so we have to
+ // copy the data vectors a bit around
+ //
+ // note that we copy vectors when
+ // looping over the patches since we
+ // have to write them one variable
+ // at a time and don't want to use
+ // more than one loop
+ //
+ // this copying of data vectors can
+ // be done while we already output
+ // the vertices, so do this on a
+ // separate thread and when wanting
+ // to write out the data, we wait
+ // for that thread to finish
+ std::vector<std::vector<double> > data_vectors (n_data_sets,
+ std::vector<double> (n_nodes));
+ Threads::ThreadManager thread_manager;
+ Threads::spawn (thread_manager,
+ Threads::encapsulate (&DataOutBase::template
+ write_gmv_reorder_data_vectors<dim,spacedim>)
+ .collect_args(patches, data_vectors));
+
+ ///////////////////////////////
+ // first make up a list of used
+ // vertices along with their
+ // coordinates
+ //
+ // note that we have to print
+ // d=1..3 dimensions
+ out << "POINTS " << n_nodes << " double" << std::endl;
+ for (unsigned int d=1; d<=3; ++d)
+ {
+ for (typename std::vector<Patch<dim,spacedim> >::const_iterator patch=patches.begin();
+ patch!=patches.end(); ++patch)
+ {
+ const unsigned int n_subdivisions = patch->n_subdivisions;
+
+ // if we have nonzero values for
+ // this coordinate
+ if (d<=spacedim)
+ {
+ switch (dim)
+ {
+ case 1:
+ {
+ for (unsigned int i=0; i<n_subdivisions+1; ++i)
+ out << ((patch->vertices[1](0) * i / n_subdivisions) +
+ (patch->vertices[0](0) * (n_subdivisions-i) / n_subdivisions))
+ << ' ';
+ break;
+ };
+
+ case 2:
+ {
+ for (unsigned int i=0; i<n_subdivisions+1; ++i)
+ for (unsigned int j=0; j<n_subdivisions+1; ++j)
+ {
+ const double x_frac = i * 1./n_subdivisions,
+ y_frac = j * 1./n_subdivisions;
+
+ // compute coordinates for
+ // this patch point
+ out << (((patch->vertices[1](d-1) * x_frac) +
+ (patch->vertices[0](d-1) * (1-x_frac))) * (1-y_frac) +
+ ((patch->vertices[2](d-1) * x_frac) +
+ (patch->vertices[3](d-1) * (1-x_frac))) * y_frac)
+ << ' ';
+ };
+ break;
+ };
+
+ case 3:
+ {
+ for (unsigned int i=0; i<n_subdivisions+1; ++i)
+ for (unsigned int j=0; j<n_subdivisions+1; ++j)
+ for (unsigned int k=0; k<n_subdivisions+1; ++k)
+ {
+ // note the broken
+ // design of hexahedra
+ // in deal.II, where
+ // first the z-component
+ // is counted up, before
+ // increasing the y-
+ // coordinate
+ const double x_frac = i * 1./n_subdivisions,
+ y_frac = k * 1./n_subdivisions,
+ z_frac = j * 1./n_subdivisions;
+
+ // compute coordinates for
+ // this patch point
+ out << ((((patch->vertices[1](d-1) * x_frac) +
+ (patch->vertices[0](d-1) * (1-x_frac))) * (1-y_frac) +
+ ((patch->vertices[2](d-1) * x_frac) +
+ (patch->vertices[3](d-1) * (1-x_frac))) * y_frac) * (1-z_frac) +
+ (((patch->vertices[5](d-1) * x_frac) +
+ (patch->vertices[4](d-1) * (1-x_frac))) * (1-y_frac) +
+ ((patch->vertices[6](d-1) * x_frac) +
+ (patch->vertices[7](d-1) * (1-x_frac))) * y_frac) * z_frac)
+ << ' ';
+ };
+
+ break;
+ };
+
+ default:
+ Assert (false, ExcNotImplemented());
+ };
+ }
+ else
+ // d>spacedim. write zeros instead
+ {
+ const unsigned int n_points
+ = static_cast<unsigned int>(pow (static_cast<double>(n_subdivisions+1), dim));
+ for (unsigned int i=0; i<n_points; ++i)
+ out << "0 ";
+ };
+ };
+ out << std::endl;
+ };
+
+ out << std::endl;
+
+ /////////////////////////////////
+ // now for the cells. note that
+ // vertices are counted from 1 onwards
+ if (true)
+ {
+ out << "CELLS " << n_cells << ' '
+ << n_cells*(GeometryInfo<dim>::vertices_per_cell+1)
+ << std::endl;
+
+
+ unsigned int first_vertex_of_patch = 0;
+
+ for (typename std::vector<Patch<dim,spacedim> >::const_iterator patch=patches.begin();
+ patch!=patches.end(); ++patch)
+ {
+ const unsigned int n_subdivisions = patch->n_subdivisions;
+
+ // write out the cells making
+ // up this patch
+ switch (dim)
+ {
+ case 1:
+ {
+ for (unsigned int i=0; i<n_subdivisions; ++i)
+ out << "2 "
+ << first_vertex_of_patch+i+1 << ' '
+ << first_vertex_of_patch+i+1+1 << std::endl;
+ break;
+ };
+
+ case 2:
+ {
+ for (unsigned int i=0; i<n_subdivisions; ++i)
+ for (unsigned int j=0; j<n_subdivisions; ++j)
+ out << "4 "
+ << first_vertex_of_patch+i*(n_subdivisions+1)+j+1 << ' '
+ << first_vertex_of_patch+(i+1)*(n_subdivisions+1)+j+1 << ' '
+ << first_vertex_of_patch+(i+1)*(n_subdivisions+1)+j+1+1 << ' '
+ << first_vertex_of_patch+i*(n_subdivisions+1)+j+1+1
+ << std::endl;
+ break;
+ };
+
+ case 3:
+ {
+ for (unsigned int i=0; i<n_subdivisions; ++i)
+ for (unsigned int j=0; j<n_subdivisions; ++j)
+ for (unsigned int k=0; k<n_subdivisions; ++k)
+ {
+ out << "8 "
+ // note: vertex indices start with 1!
+ << first_vertex_of_patch+(i*(n_subdivisions+1)+j )*(n_subdivisions+1)+k +1 << ' '
+ << first_vertex_of_patch+((i+1)*(n_subdivisions+1)+j )*(n_subdivisions+1)+k +1 << ' '
+ << first_vertex_of_patch+((i+1)*(n_subdivisions+1)+j+1)*(n_subdivisions+1)+k +1 << ' '
+ << first_vertex_of_patch+(i*(n_subdivisions+1)+j+1 )*(n_subdivisions+1)+k +1 << ' '
+ << first_vertex_of_patch+(i*(n_subdivisions+1)+j )*(n_subdivisions+1)+k+1+1 << ' '
+ << first_vertex_of_patch+((i+1)*(n_subdivisions+1)+j )*(n_subdivisions+1)+k+1+1 << ' '
+ << first_vertex_of_patch+((i+1)*(n_subdivisions+1)+j+1)*(n_subdivisions+1)+k+1+1 << ' '
+ << first_vertex_of_patch+(i*(n_subdivisions+1)+j+1 )*(n_subdivisions+1)+k+1+1 << ' '
+ << std::endl;
+ };
+ break;
+ };
+
+ default:
+ Assert (false, ExcNotImplemented());
+ };
+
+
+ // finally update the number
+ // of the first vertex of this patch
+ switch (dim)
+ {
+ case 1:
+ first_vertex_of_patch += n_subdivisions+1;
+ break;
+ case 2:
+ first_vertex_of_patch += (n_subdivisions+1) *
+ (n_subdivisions+1);
+ break;
+ case 3:
+ first_vertex_of_patch += (n_subdivisions+1) *
+ (n_subdivisions+1) *
+ (n_subdivisions+1);
+ break;
+ default:
+ Assert (false, ExcNotImplemented());
+ };
+ };
+ out << std::endl;
+
+ // next output the types of the
+ // cells. since all cells are
+ // the same, this is simple
+ for (unsigned int i=0; i<n_cells; ++i)
+ switch (dim)
+ {
+ case 1:
+ out << "3\n"; // represents VTK_LINE
+ break;
+ case 2:
+ out << "9\n"; // represents VTK_QUAD
+ break;
+ case 3:
+ out << "12\n"; // represents VTK_HEXAHEDRON
+ break;
+ default:
+ Assert (false, ExcNotImplemented());
+ };
+ };
+
+ ///////////////////////////////////////
+ // data output.
+
+ // now write the data vectors to
+ // @p{out} first make sure that all
+ // data is in place
+ thread_manager.wait ();
+
+ // then write data.
+ // the '1' means: node data (as opposed
+ // to cell data, which we do not
+ // support explicitely here)
+ for (unsigned int data_set=0; data_set<n_data_sets; ++data_set)
+ {
+ out << "POINTDATA " << n_nodes
+ << std::endl
+ << "SCALARS "
+ << data_names[data_set]
+ << " double 1"
+ << std::endl
+ << "LOOKUP_TABLE default"
+ << std::endl;
+ copy(data_vectors[data_set].begin(),
+ data_vectors[data_set].end(),
+ std::ostream_iterator<double>(out, " "));
+ out << std::endl
+ << std::endl;
+ };
+
+ // assert the stream is still ok
+ AssertThrow (out, ExcIO());
+};
+
+
+
+
template <int dim, int spacedim>
void
DataOutBase::write_gmv_reorder_data_vectors (const typename std::vector<Patch<dim,spacedim> > &patches,
+template <int dim, int spacedim>
+void DataOutInterface<dim,spacedim>::write_vtk (std::ostream &out) const
+{
+ DataOutBase::write_vtk (get_patches(), get_dataset_names(),
+ vtk_flags, out);
+};
+
+
+
template <int dim, int spacedim>
void DataOutInterface<dim,spacedim>::write (std::ostream &out,
OutputFormat output_format) const
write_gmv (out);
break;
+ case vtk:
+ write_vtk (out);
+ break;
+
default:
Assert (false, ExcNotImplemented());
};
+template <int dim, int spacedim>
+void DataOutInterface<dim,spacedim>::set_flags (const VtkFlags &flags)
+{
+ vtk_flags = flags;
+};
+
+
+
template <int dim, int spacedim>
std::string DataOutInterface<dim,spacedim>::default_suffix (OutputFormat output_format) const
{
case gmv:
return ".gmv";
+ case vtk:
+ return ".vtk";
+
default:
Assert (false, ExcNotImplemented());
return "";
if (format_name == "gmv")
return gmv;
+ if (format_name == "vtk")
+ return vtk;
+
AssertThrow (false, ExcInvalidState ());
// return something invalid
template <int dim, int spacedim>
std::string DataOutInterface<dim,spacedim>::get_output_format_names ()
{
- return "ucd|gnuplot|povray|eps|gmv";
+ return "ucd|gnuplot|povray|eps|gmv|vtk";
}
prm.enter_subsection ("Gmv output parameters");
GmvFlags::declare_parameters (prm);
prm.leave_subsection ();
+
+ prm.enter_subsection ("Vtk output parameters");
+ VtkFlags::declare_parameters (prm);
+ prm.leave_subsection ();
}
prm.enter_subsection ("Gmv output parameters");
gmv_flags.parse_parameters (prm);
prm.leave_subsection ();
+
+ prm.enter_subsection ("Vtk output parameters");
+ vtk_flags.parse_parameters (prm);
+ prm.leave_subsection ();
}
MemoryConsumption::memory_consumption (gnuplot_flags) +
MemoryConsumption::memory_consumption (povray_flags) +
MemoryConsumption::memory_consumption (eps_flags) +
- MemoryConsumption::memory_consumption (gmv_flags));
+ MemoryConsumption::memory_consumption (gmv_flags) +
+ MemoryConsumption::memory_consumption (vtk_flags));
};