]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add .visit time information 3012/head
authorRene Gassmoeller <rene.gassmoeller@mailbox.org>
Wed, 24 Aug 2016 23:11:34 +0000 (17:11 -0600)
committerRene Gassmoeller <rene.gassmoeller@mailbox.org>
Fri, 26 Aug 2016 00:14:01 +0000 (18:14 -0600)
doc/news/changes.h
include/deal.II/base/data_out_base.h
source/base/data_out_base.cc
tests/numerics/data_out_10.cc [new file with mode: 0644]
tests/numerics/data_out_10.output [new file with mode: 0644]

index d255170027a37104cdd0c3406e207595964d6193..690883db3040130b9f329180f9c091a1ad3cc1cb 100644 (file)
@@ -38,6 +38,13 @@ inconvenience this causes.
 </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
index 0157e5e41f898b37b8fe48c9e2496d09acab6025..78eb6e661912cc9995f74252fd093cce23a3dbfc 100644 (file)
@@ -2225,6 +2225,42 @@ public:
   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> > > &times_and_piece_names) const;
+
   /**
    * Obtain data through get_patches() and write it to <tt>out</tt> in SVG
    * format. See DataOutBase::write_svg.
index 3fbe0e9837d05fcfd70be015ef130e84c7c983f9..43d0829517963c975443a5e9e2881469671faa3c 100644 (file)
@@ -6384,6 +6384,37 @@ DataOutInterface<dim,spacedim>::write_visit_record (std::ostream &out,
 
 
 
+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> > > &times_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
diff --git a/tests/numerics/data_out_10.cc b/tests/numerics/data_out_10.cc
new file mode 100644 (file)
index 0000000..b45bb69
--- /dev/null
@@ -0,0 +1,70 @@
+// ---------------------------------------------------------------------
+//
+// 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();
+
+}
+
diff --git a/tests/numerics/data_out_10.output b/tests/numerics/data_out_10.output
new file mode 100644 (file)
index 0000000..25888a4
--- /dev/null
@@ -0,0 +1,12 @@
+
+!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

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.