* there for more information.
*
*
+ * @subsection DataOutBaseD2 deal.II intermediate format
+ *
+ * In addition to all the other formats, this class can also write
+ * data in the deal.II format. This is not a format that is read 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. The reason why we offer to write out this intermediate
+ * format is that it can be read back into a deal.II program, 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 <tt>.d2</tt>.
+ *
+ *
* @section DataOutBaseOP Output parameters
*
* All functions take a parameter which is a structure of type <tt>XFlags</tt>, where
unsigned int memory_consumption () const;
/**
- * Value for no neighbor.
+ * Value to be used if this
+ * patch has no neighbor on
+ * one side.
*/
static const unsigned int no_neighbor = deal_II_numbers::invalid_unsigned_int;
/** @addtogroup Exceptions
unsigned int memory_consumption () const;
};
+
+
+ /**
+ * Flags describing the details
+ * of output in deal.II
+ * intermediate format. At
+ * present no flags are
+ * implemented.
+ */
+ struct Deal_II_IntermediateFlags
+ {
+ 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
+ * <tt>private</tt> 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
+ * <tt>declare_parameters</tt> 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 <tt>std::map</tt> 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;
+ };
+
/**
* Write the given list of patches
* to the output stream in OpenDX
const VtkFlags &flags,
std::ostream &out);
+ /**
+ * Write the given list of
+ * patches to the output stream
+ * in deal.II intermediate
+ * format. See the general
+ * documentation for more
+ * information on the parameters.
+ */
+ template <int dim, int spacedim>
+ static void write_deal_II_intermediate (const std::vector<Patch<dim,spacedim> > &patches,
+ const std::vector<std::string> &data_names,
+ const Deal_II_IntermediateFlags &flags,
+ std::ostream &out);
+
/**
* Determine an estimate for
* the memory consumption (in
/**
* Output in VTK format.
*/
- vtk };
+ vtk,
+ /**
+ * Output in deal.II
+ * intermediate format.
+ */
+ deal_II_intermediate
+ };
/**
* Destructor. Does nothing, but is
*/
void write_vtk (std::ostream &out) const;
+ /**
+ * Obtain data through the
+ * <tt>get_patches</tt> function
+ * and write it to <tt>out</tt>
+ * in deal.II intermediate
+ * format.
+ */
+ void write_deal_II_intermediate (std::ostream &out) const;
+
/**
* Write data and grid to <tt>out</tt>
* according to the given data
/**
* Set the flags to be used for
- * output in GMV format.
+ * output in VTK format.
*/
void set_flags (const VtkFlags &vtk_flags);
+
+ /**
+ * Set the flags to be used for
+ * output in VTK format.
+ */
+ void set_flags (const Deal_II_IntermediateFlags &deal_II_intermediate_flags);
/**
* <li> <tt>eps</tt>: <tt>.eps</tt>
* <li> <tt>gmv</tt>: <tt>.gmv</tt>
* <li> <tt>vtk</tt>: <tt>.vtk</tt>.
+ * <li> <tt>deal_II_intermediate</tt>: <tt>.d2</tt>.
* </ul>
*
* If this function is called
* function.
*/
VtkFlags vtk_flags;
+
+ /**
+ * Flags to be used upon output
+ * of vtk data in one space
+ * dimension. Can be changed by
+ * using the <tt>set_flags</tt>
+ * function.
+ */
+ Deal_II_IntermediateFlags deal_II_intermediate_flags;
};
}
+/* -------------------- template functions ------------------- */
+
+/**
+ * Output operator for an object of type
+ * <tt>DataOutBase::Patch</tt>. 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.
+ *
+ * @author Wolfgang Bangerth, 2005
+ */
+template <int dim, int spacedim>
+std::ostream &
+operator << (std::ostream &out,
+ const DataOutBase::Patch<dim,spacedim> &patch);
+
+
#endif
+void DataOutBase::Deal_II_IntermediateFlags::declare_parameters (ParameterHandler &/*prm*/)
+{}
+
+
+
+void DataOutBase::Deal_II_IntermediateFlags::parse_parameters (ParameterHandler &/*prm*/)
+{}
+
+
+unsigned int
+DataOutBase::Deal_II_IntermediateFlags::memory_consumption () const
+{
+ // only simple data elements, so
+ // use sizeof operator
+ return sizeof (*this);
+}
+
+
+
unsigned int DataOutBase::memory_consumption ()
{
return 0;
+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 Deal_II_IntermediateFlags &/*flags*/,
+ std::ostream &out)
+{
+ out << "[deal.II intermediate format graphics data]" << std::endl
+ << "[written by deal.II version "
+ << DEAL_II_MAJOR << '.' << DEAL_II_MINOR << "]" << std::endl;
+
+ out << data_names.size() << std::endl;
+ for (unsigned int i=0; i<data_names.size(); ++i)
+ out << data_names[i] << std::endl;
+
+ out << patches.size() << std::endl;
+ for (unsigned int i=0; i<patches.size(); ++i)
+ out << patches[i] << std::endl;
+
+ out << std::endl;
+}
+
+
template <int dim, int spacedim>
void
+template <int dim, int spacedim>
+void DataOutInterface<dim,spacedim>::
+write_deal_II_intermediate (std::ostream &out) const
+{
+ DataOutBase::write_deal_II_intermediate (get_patches(), get_dataset_names(),
+ deal_II_intermediate_flags, out);
+}
+
+
+
template <int dim, int spacedim>
void
DataOutInterface<dim,spacedim>::write (std::ostream &out,
case vtk:
write_vtk (out);
break;
+
+ case deal_II_intermediate:
+ write_deal_II_intermediate (out);
+ break;
default:
Assert (false, ExcNotImplemented());
+template <int dim, int spacedim>
+void
+DataOutInterface<dim,spacedim>::set_flags (const Deal_II_IntermediateFlags &flags)
+{
+ deal_II_intermediate_flags = flags;
+}
+
+
+
template <int dim, int spacedim>
std::string
DataOutInterface<dim,spacedim>::
case vtk:
return ".vtk";
+ case deal_II_intermediate:
+ return ".d2";
+
default:
Assert (false, ExcNotImplemented());
return "";
if (format_name == "vtk")
return vtk;
+ if (format_name == "deal.II intermediate")
+ return deal_II_intermediate;
+
AssertThrow (false, ExcInvalidState ());
// return something invalid
std::string
DataOutInterface<dim,spacedim>::get_output_format_names ()
{
- return "dx|ucd|gnuplot|povray|eps|gmv|tecplot|vtk";
+ return "dx|ucd|gnuplot|povray|eps|gmv|tecplot|vtk|deal.II intermediate";
}
prm.enter_subsection ("Vtk output parameters");
VtkFlags::declare_parameters (prm);
prm.leave_subsection ();
+
+
+ prm.enter_subsection ("deal.II intermediate output parameters");
+ Deal_II_IntermediateFlags::declare_parameters (prm);
+ prm.leave_subsection ();
}
prm.enter_subsection ("Vtk output parameters");
vtk_flags.parse_parameters (prm);
prm.leave_subsection ();
+
+ prm.enter_subsection ("deal.II intermediate output parameters");
+ deal_II_intermediate_flags.parse_parameters (prm);
+ prm.leave_subsection ();
}
MemoryConsumption::memory_consumption (eps_flags) +
MemoryConsumption::memory_consumption (gmv_flags) +
MemoryConsumption::memory_consumption (tecplot_flags) +
- MemoryConsumption::memory_consumption (vtk_flags));
+ MemoryConsumption::memory_consumption (vtk_flags) +
+ MemoryConsumption::memory_consumption (deal_II_intermediate_flags));
+}
+
+
+
+template <int dim, int spacedim>
+std::ostream &
+operator << (std::ostream &out,
+ const DataOutBase::Patch<dim,spacedim> &patch)
+{
+ // write a header line
+ out << "[deal.II intermediate Patch<" << dim << ',' << spacedim << ">]"
+ << std::endl;
+
+ // then write all the data that is
+ // in this patch
+ for (unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i)
+ out << patch.vertices[i] << " ";
+ out << std::endl;
+
+ for (unsigned int i=0; i<GeometryInfo<dim>::faces_per_cell; ++i)
+ out << patch.neighbors[i] << " ";
+ out << std::endl;
+
+ out << patch.patch_index << ' ' << patch.n_subdivisions
+ << std::endl;
+
+ out << patch.data.n_rows() << ' ' << patch.data.n_cols() << std::endl;
+ for (unsigned int i=0; i<patch.data.n_rows(); ++i)
+ for (unsigned int j=0; j<patch.data.n_cols(); ++j)
+ out << patch.data[i][j] << ' ';
+ out << std::endl;
+ out << std::endl;
+
+ return out;
}
template class DataOutInterface<3,4>;
template class DataOutBase::Patch<3,4>;
+// output operators
+template
+std::ostream &
+operator << (std::ostream &out,
+ const DataOutBase::Patch<1,1> &patch);
+template
+std::ostream &
+operator << (std::ostream &out,
+ const DataOutBase::Patch<2,2> &patch);
+template
+std::ostream &
+operator << (std::ostream &out,
+ const DataOutBase::Patch<3,3> &patch);
+template
+std::ostream &
+operator << (std::ostream &out,
+ const DataOutBase::Patch<4,4> &patch);
+template
+std::ostream &
+operator << (std::ostream &out,
+ const DataOutBase::Patch<1,2> &patch);
+template
+std::ostream &
+operator << (std::ostream &out,
+ const DataOutBase::Patch<2,3> &patch);
+template
+std::ostream &
+operator << (std::ostream &out,
+ const DataOutBase::Patch<3,4> &patch);
+
+