]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Patch by Marco Engelhard: Add DataOutInterface::write_pvd_record.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 30 May 2012 12:53:31 +0000 (12:53 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 30 May 2012 12:53:31 +0000 (12:53 +0000)
git-svn-id: https://svn.dealii.org/trunk@25569 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/base/data_out_base.h
deal.II/source/base/data_out_base.cc
tests/base/data_out_base_pvd.cc [new file with mode: 0644]
tests/base/data_out_base_pvd/cmp/generic [new file with mode: 0644]

index 68f51af9530d64a04c30a387002a9629848cb430..be820c9c3b7b1bd38e6761d53b6f5b98bee2d6bb 100644 (file)
@@ -48,12 +48,6 @@ used for boundary indicators.
 
 <ol>
 
-<li> New: DoFTools::make_periodicity_constraints implemented which
-inserts algebraic constraints due to periodic boundary conditions
-into a ConstraintMatrix.
-<br>
-(Matthias Maier, 2012/05/22)
-
 <li>
 New: step-48 demonstrates the solution of a nonlinear wave equation
 with an explicit time stepping method. The usage of Gauss-Lobatto
@@ -81,7 +75,7 @@ implemented using the FEEvaluation class.
 (Katharina Kormann and Martin Kronbichler, 2012/05/05)
 
 <li>
-New: step-44 demonstrates one approach to modelling large
+New: step-44 demonstrates one approach to modeling large
 deformations of nearly-incompressible elastic solids. The
 elastic response is governed by a non-linear, hyperelastic
 free-energy function. The geometrical response is also
@@ -207,11 +201,23 @@ enabled due to a missing include file in file
 <h3>Specific improvements</h3>
 
 <ol>
+<li> New: The function DataOutInterface::write_pvd_record can be used
+to provide Paraview with metadata that describes which time in a
+simulation a particular output file corresponds to.
+<br>
+(Marco Engelhard 2012/05/30)
+
 <li> Fixed: Bug in 3d with hanging nodes in GridTools::find_cells_adjacent_to_vertex()
 that caused find_active_cell_around_point() to fail in those cases.
 <br>
 (Timo Heister, Wolfgang Bangerth 2012/05/30)
 
+<li> New: DoFTools::make_periodicity_constraints implemented which
+inserts algebraic constraints due to periodic boundary conditions
+into a ConstraintMatrix.
+<br>
+(Matthias Maier, 2012/05/22)
+
 <li> New: The GridIn::read_unv function can now read meshes generated
 by the Salome framework, see http://www.salome-platform.org/ .
 <br>
index 560e09e55e3e9e1f570ff4be7364ae4a97c0057d..f0fc486f7402aecf0e49f0f41120e8cd41861136 100644 (file)
@@ -2348,6 +2348,61 @@ class DataOutInterface : private DataOutBase
     void write_pvtu_record (std::ostream &out,
                             const std::vector<std::string> &piece_names) const;
 
+    /**
+     * 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 <code>.pvd</code>
+     * 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 (<code>time</code> and
+     * <code>time_step</code> are member variables of the class with types
+     * <code>double</code> and <code>unsigned int</code>, respectively;
+     * the variable <code>times_and_names</code> is of type
+     * <code>std::vector@<std::pair@<double,std::string@> @></code>):
+     *
+     * @code
+     * 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_number, 3) +
+     *                                 ".vtu";
+     *    std::ofstream output (filename.c_str());
+     *    data_out.write_vtu (output);
+     *
+     *    times_and_names.push_back (std::pair<double,std::string> (time, filename));
+     *    std::ofstream pvd_output ("solution.pvd");
+     *    data_out.write_pvd_record (pvd_output, times_and_names);
+     *  }
+     * @endcode
+     *
+     * @note See DataOutBase::write_vtu or
+     * DataOutBase::write_pvtu_record for
+     * writing solutions at each timestep.
+     *
+     * @note 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 <code>.pvtu</code> file that references multiple
+     * parts of a parallel computation.
+     *
+     * @author Marco Engelhard, 2012
+     */
+    void write_pvd_record (std::ostream &out,
+                           const std::vector<std::pair<double,std::string> >  &times_and_names) const;
+
                                      /**
                                       * This function is the exact
                                       * equivalent of the
index 557c59aac0325ea530a5e0e01614214983af5e1e..98380bf5689187a60d31438e3f7d69d23fff4219 100644 (file)
@@ -5389,6 +5389,44 @@ void DataOutInterface<dim,spacedim>::write_vtu_in_parallel (const char* filename
 }
 
 
+template <int dim, int spacedim>
+void
+DataOutInterface<dim,spacedim>::
+write_pvd_record (std::ostream &out,
+                  const std::vector<std::pair<double,std::string> >  &times_and_names) const
+{
+  AssertThrow (out, ExcIO());
+
+  out << "<?xml version=\"1.0\"?>\n";
+
+  std::time_t  time1= std::time (0);
+  std::tm     *time = std::localtime(&time1);
+  out << "<!--\n";
+  out << "#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
+      << "\n-->\n";
+
+  out << "<VTKFile type=\"Collection\" version=\"0.1\" ByteOrder=\"LittleEndian\">\n";
+  out << "  <Collection>\n";
+
+  for(unsigned int i=0; i<times_and_names.size(); ++i)
+    out << "    <DataSet timestep=\"" << times_and_names[i].first
+        << "\" group=\"\" part=\"0\" file=\"" << times_and_names[i].second
+        << "\"/>\n";
+
+  out << "  </Collection>\n";
+  out << "</VTKFile>\n";
+
+  out.flush();
+
+  AssertThrow (out, ExcIO());
+}
+
 
 template <int dim, int spacedim>
 void
diff --git a/tests/base/data_out_base_pvd.cc b/tests/base/data_out_base_pvd.cc
new file mode 100644 (file)
index 0000000..8af1b77
--- /dev/null
@@ -0,0 +1,71 @@
+//-----------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2006, 2007, 2010, 2012 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//-----------------------------------------------------------------------------
+
+// write the pvd master record for parallel visualization through the
+// vtu file format
+
+#include "../tests.h"
+#include <deal.II/base/data_out_base.h>
+#include <deal.II/base/logstream.h>
+
+#include <vector>
+#include <iomanip>
+#include <fstream>
+#include <string>
+#include <stdio.h>
+
+#include "patches.h"
+
+
+
+std::vector<DataOutBase::Patch<2,2> > patches;
+
+class DataOutX : public DataOutInterface<2,2>
+{
+    virtual
+    const std::vector< ::DataOutBase::Patch<2,2> > &
+    get_patches () const
+      {
+       return patches;
+      }
+
+    virtual
+    std::vector<std::string>
+    get_dataset_names () const
+      {
+       return std::vector<std::string>();
+      }
+};
+
+
+template <int dim, int spacedim>
+void check(std::ostream& out)
+{
+  std::vector<std::pair<double,std::string> > names(5);
+  names[0] = std::make_pair(0,"x1");
+  names[1] = std::make_pair(1,"x2");
+  names[2] = std::make_pair(1e1,"x3");
+  names[3] = std::make_pair(3.141,"d");
+  names[4] = std::make_pair(42e19,"i");
+
+  DataOutX x;
+  x.write_pvd_record (out, names);
+}
+
+
+
+int main()
+{
+  std::ofstream logfile("data_out_base_pvd/output");
+  check<2,2>(logfile);
+}
diff --git a/tests/base/data_out_base_pvd/cmp/generic b/tests/base/data_out_base_pvd/cmp/generic
new file mode 100644 (file)
index 0000000..8524c77
--- /dev/null
@@ -0,0 +1,13 @@
+<?xml version="1.0"?>
+<!--
+#This file was generated 
+-->
+<VTKFile type="Collection" version="0.1" ByteOrder="LittleEndian">
+  <Collection>
+    <DataSet timestep="0" group="" part="0" file="x1"/>
+    <DataSet timestep="1" group="" part="0" file="x2"/>
+    <DataSet timestep="10" group="" part="0" file="x3"/>
+    <DataSet timestep="3.141" group="" part="0" file="d"/>
+    <DataSet timestep="4.2e+20" group="" part="0" file="i"/>
+  </Collection>
+</VTKFile>

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.