</p>
<ol>
+ <li> New: There is now the possibility to store information about the
+ time of an output time step within the .visit file created by
+ the DataOutInterface<dim,spacedim>::write_visit_record function.
+ <br>
+ (Rene Gassmoeller, Juliane Dannberg, 2016/08/24)
+ </li>
+
<li> New: deal.II now requires at least BOOST version 1.56, rather than the
previous minimal version of 1.54. This is because 1.54 does not support
serializing objects of type std::unique_ptr if C++11 is used, but we now
void write_visit_record (std::ostream &out,
const std::vector<std::vector<std::string> > &piece_names) const;
+ /**
+ * 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:
+ * @code
+ * DataOut<dim> data_out;
+ *
+ * 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.push_back("subdomain_01.time_step_0.vtk");
+ * times_and_piece_names[0].second.push_back("subdomain_02.time_step_0.vtk");
+ *
+ * times_and_piece_names[1].first = 0.5;
+ * times_and_piece_names[1].second.push_back("subdomain_01.time_step_1.vtk");
+ * times_and_piece_names[1].second.push_back("subdomain_02.time_step_1.vtk");
+ *
+ * times_and_piece_names[2].first = 1.0;
+ * times_and_piece_names[2].second.push_back("subdomain_01.time_step_2.vtk");
+ * times_and_piece_names[2].second.push_back("subdomain_02.time_step_2.vtk");
+ *
+ * std::ofstream visit_output ("master_file.visit");
+ *
+ * data_out.write_visit_record(visit_output, times_and_piece_names);
+ * @endcode
+ *
+ * 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
+ */
+ void write_visit_record (std::ostream &out,
+ const std::vector<std::pair<double,std::vector<std::string> > > ×_and_piece_names) const;
+
/**
* Obtain data through get_patches() and write it to <tt>out</tt> in SVG
* format. See DataOutBase::write_svg.
+template <int dim, int spacedim>
+void
+DataOutInterface<dim,spacedim>::write_visit_record (std::ostream &out,
+ const std::vector<std::pair<double,std::vector<std::string> > > ×_and_piece_names) const
+{
+ AssertThrow (out, ExcIO());
+
+ if (times_and_piece_names.size() == 0)
+ return;
+
+ const double nblocks = times_and_piece_names[0].second.size();
+ Assert(nblocks > 0, ExcMessage("time_and_piece_names should contain nonempty vectors of filenames for every timestep.") )
+
+ for (std::vector<std::pair<double,std::vector<std::string> > >::const_iterator domain = times_and_piece_names.begin();
+ domain != times_and_piece_names.end(); ++domain)
+ out << "!TIME " << domain->first << '\n';
+
+ out << "!NBLOCKS " << nblocks << '\n';
+ for (std::vector<std::pair<double,std::vector<std::string> > >::const_iterator domain = times_and_piece_names.begin();
+ domain != times_and_piece_names.end(); ++domain)
+ {
+ Assert(domain->second.size() == nblocks, ExcMessage("piece_names should be a vector of equal sized vectors.") )
+ for (std::vector<std::string>::const_iterator subdomain = domain->second.begin(); subdomain != domain->second.end(); ++subdomain)
+ out << *subdomain << '\n';
+ }
+
+ out << std::flush;
+}
+
+
+
template <int dim, int spacedim>
void DataOutInterface<dim,spacedim>::
write_deal_II_intermediate (std::ostream &out) const
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+// Check that adding time information to a .visit file works as intended
+
+
+#include "../tests.h"
+#include <deal.II/lac/sparsity_pattern.h>
+#include <deal.II/numerics/data_out.h>
+
+#include <deal.II/base/logstream.h>
+#include <deal.II/lac/vector.h>
+#include <deal.II/lac/block_vector.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_tools.h>
+#include <deal.II/fe/fe_dgq.h>
+
+
+
+
+void test()
+{
+ DataOut<2> data_out;
+
+ 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.push_back("subdomain-01.time_step_0.vtk");
+ times_and_piece_names[0].second.push_back("subdomain-02.time_step_0.vtk");
+
+ times_and_piece_names[1].first = 0.5;
+ times_and_piece_names[1].second.push_back("subdomain-01.time_step_1.vtk");
+ times_and_piece_names[1].second.push_back("subdomain-02.time_step_1.vtk");
+
+ times_and_piece_names[2].first = 1.0;
+ times_and_piece_names[2].second.push_back("subdomain-01.time_step_2.vtk");
+ times_and_piece_names[2].second.push_back("subdomain-02.time_step_2.vtk");
+ data_out.write_visit_record(deallog.get_file_stream(), times_and_piece_names);
+
+ deallog << "OK" << std::endl;
+
+}
+
+
+int main()
+{
+ deal_II_exceptions::disable_abort_on_exception();
+ initlog();
+ test();
+
+}
+
--- /dev/null
+
+!TIME 0.00000
+!TIME 0.500000
+!TIME 1.00000
+!NBLOCKS 2.00000
+subdomain-01.time_step_0.vtk
+subdomain-02.time_step_0.vtk
+subdomain-01.time_step_1.vtk
+subdomain-02.time_step_1.vtk
+subdomain-01.time_step_2.vtk
+subdomain-02.time_step_2.vtk
+DEAL::OK