]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
New: DataOut::write_vtu_in_parallel(). This routine uses MPI I/O to write a big vtu...
authorheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 1 Mar 2012 00:10:18 +0000 (00:10 +0000)
committerheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 1 Mar 2012 00:10:18 +0000 (00:10 +0000)
git-svn-id: https://svn.dealii.org/trunk@25195 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

index d30b5b3b8b902764bdcc5c06042a3da9bff5985c..177f59bbe271fbf23d3128616b02caf40cf3bd12 100644 (file)
@@ -145,6 +145,11 @@ enabled due to a missing include file in file
 <h3>Specific improvements</h3>
 
 <ol>
+<li> New: DataOut::write_vtu_in_parallel(). This routine uses MPI I/O to write
+a big vtu file in parallel.
+<br>
+(Timo Heister, 2012/02/29)
+
 <li> Fixed: parallel::distributed::Triangulation::clear() forgot
 to update the number cache of this class, leading to wrong results
 when calling functions like
index 887a23dd431c047c27dd3669188dbe5507b630c0..077a3e5622ebe94d4fa57a19923986fa5f3329ba 100644 (file)
@@ -22,6 +22,8 @@
 #include <vector>
 #include <string>
 
+#include <deal.II/base/mpi.h>
+
 // Only include the Tecplot API header if the appropriate files
 // were detected by configure
 #ifdef DEAL_II_HAVE_TECPLOT
@@ -1723,6 +1725,35 @@ class DataOutBase
                           const VtkFlags                          &flags,
                           std::ostream                            &out);
 
+/**
+ * This writes the header for the xml based vtu file format. This
+ * routine is used internally together with
+ * DataOutInterface::write_vtu_footer() and DataOutInterface::write_vtu_main()
+ * by DataOutBase::write_vtu().
+ */
+    static void write_vtu_header (std::ostream &out);
+
+/**
+ * This writes the footer for the xml based vtu file format. This
+ * routine is used internally together with
+ * DataOutInterface::write_vtu_header() and DataOutInterface::write_vtu_main()
+ * by DataOutBase::write_vtu().
+ */
+    static void write_vtu_footer (std::ostream &out);
+
+/**
+ * This writes the main part for the xml based vtu file format. This
+ * routine is used internally together with
+ * DataOutInterface::write_vtu_header() and DataOutInterface::write_vtu_footer()
+ * by DataOutBase::write_vtu().
+ */
+    template <int dim, int spacedim>
+    static void write_vtu_main (const std::vector<Patch<dim,spacedim> > &patches,
+                         const std::vector<std::string>          &data_names,
+                         const std::vector<std_cxx1x::tuple<unsigned int, unsigned int, std::string> > &vector_data_ranges,
+                         const VtkFlags                          &flags,
+                         std::ostream                            &out);
+
 /**
  * Write the given list of patches to the output stream in deal.II
  * intermediate format. This is not a format understood by any other
@@ -1995,6 +2026,8 @@ class DataOutBase
     write_gmv_reorder_data_vectors (const std::vector<Patch<dim,spacedim> > &patches,
                                    Table<2,double>                         &data_vectors);
 
+
+
 };
 
 
@@ -2253,6 +2286,21 @@ class DataOutInterface : private DataOutBase
                                      */
     void write_vtu (std::ostream &out) const;
 
+                                     /**
+                                      * Collective MPI call to write the
+                                     * solution from all participating nodes
+                                     * (those in the given communicator) to a
+                                     * single compressed .vtu file on a
+                                     * shared file system.  The communicator
+                                     * can be a sub communicator of the one
+                                     * used by the computation.  This routine
+                                     * uses MPI I/O to achieve high
+                                     * performance on parallel filesystems.
+                                     * Also see
+                                     * DataOutInterface::write_vtu().
+                                      */
+    void write_vtu_in_parallel (const char* filename, MPI_Comm comm) const;
+
                                     /**
                                      * Some visualization programs, such as
                                      * ParaView, can read several separate
index 70290c412f682133568b0b7850b97a74d0e8c0af..ebbf91020feb0bde9e74cf053881c23a4048cd71 100644 (file)
@@ -34,6 +34,7 @@
 #include <deal.II/base/thread_management.h>
 #include <deal.II/base/memory_consumption.h>
 #include <deal.II/base/std_cxx1x/shared_ptr.h>
+#include <deal.II/base/mpi.h>
 
 #include <cstring>
 #include <algorithm>
@@ -41,8 +42,8 @@
 #include <ctime>
 #include <cmath>
 #include <set>
-
 #include <sstream>
+#include <fstream>
 
 #ifdef HAVE_LIBZ
 #  include <zlib.h>
@@ -4710,13 +4711,66 @@ DataOutBase::write_vtk (const std::vector<Patch<dim,spacedim> > &patches,
 }
 
 
+void DataOutBase::write_vtu_header (std::ostream &out)
+{
+  AssertThrow (out, ExcIO());
+  std::time_t  time1= std::time (0);
+  std::tm     *time = std::localtime(&time1);
+  out << "<?xml version=\"1.0\" ?> \n";
+  out << "<!-- \n";
+  out << "# vtk DataFile Version 3.0"
+      << '\n'
+      << "#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=\"UnstructuredGrid\" version=\"0.1\"";
+#ifdef HAVE_LIBZ
+  out << " compressor=\"vtkZLibDataCompressor\"";
+#endif
+#ifdef DEAL_II_WORDS_BIGENDIAN
+  out << " byte_order=\"BigEndian\"";
+#else
+  out << " byte_order=\"LittleEndian\"";
+#endif
+  out << ">";
+  out << '\n';
+  out << "<UnstructuredGrid>";
+  out << '\n';
+}
+
+void DataOutBase::write_vtu_footer (std::ostream &out)
+{
+  AssertThrow (out, ExcIO());
+  out << " </UnstructuredGrid>\n";
+  out << "</VTKFile>\n";
+}
+
 template <int dim, int spacedim>
 void
 DataOutBase::write_vtu (const std::vector<Patch<dim,spacedim> > &patches,
-                       const std::vector<std::string>          &data_names,
-                       const std::vector<std_cxx1x::tuple<unsigned int, unsigned int, std::string> > &vector_data_ranges,
-                       const VtkFlags                          &flags,
-                       std::ostream                            &out)
+                        const std::vector<std::string>          &data_names,
+                        const std::vector<std_cxx1x::tuple<unsigned int, unsigned int, std::string> > &vector_data_ranges,
+                        const VtkFlags                          &flags,
+                        std::ostream                            &out)
+{
+  write_vtu_header(out);
+  write_vtu_main (patches, data_names, vector_data_ranges, flags, out);
+  write_vtu_footer(out);
+}
+
+
+template <int dim, int spacedim>
+void DataOutBase::write_vtu_main (const std::vector<Patch<dim,spacedim> > &patches,
+                                  const std::vector<std::string>          &data_names,
+                                  const std::vector<std_cxx1x::tuple<unsigned int, unsigned int, std::string> > &vector_data_ranges,
+                                  const VtkFlags                          &flags,
+                                  std::ostream                            &out)
 {
   AssertThrow (out, ExcIO());
 
@@ -4756,40 +4810,6 @@ DataOutBase::write_vtu (const std::vector<Patch<dim,spacedim> > &patches,
                                n_data_sets,
                                patches[0].data.n_rows()));
 
-                                  ///////////////////////
-                                  // preamble
-  if (true)
-    {
-      std::time_t  time1= std::time (0);
-      std::tm     *time = std::localtime(&time1);
-      out << "<?xml version=\"1.0\" ?> \n";
-      out << "<!-- \n";
-      out << "# vtk DataFile Version 3.0"
-         << '\n'
-         << "#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=\"UnstructuredGrid\" version=\"0.1\"";
-#ifdef HAVE_LIBZ
-      out << " compressor=\"vtkZLibDataCompressor\"";
-#endif
-#ifdef DEAL_II_WORDS_BIGENDIAN
-      out << " byte_order=\"BigEndian\"";
-#else
-      out << " byte_order=\"LittleEndian\"";
-#endif
-      out << ">";
-      out << '\n';
-      out << "<UnstructuredGrid>";
-      out << '\n';
-    }
-
 
 #ifdef HAVE_LIBZ
   const char *ascii_or_binary = "binary";
@@ -5019,8 +5039,7 @@ DataOutBase::write_vtu (const std::vector<Patch<dim,spacedim> > &patches,
 
                                   // Finish up writing a valid XML file
   out << " </Piece>\n";
-  out << " </UnstructuredGrid>\n";
-  out << "</VTKFile>\n";
+
                                   // make sure everything now gets to
                                   // disk
   out.flush ();
@@ -5254,6 +5273,66 @@ void DataOutInterface<dim,spacedim>::write_vtu (std::ostream &out) const
                          vtk_flags, out);
 }
 
+template <int dim, int spacedim>
+void DataOutInterface<dim,spacedim>::write_vtu_in_parallel (const char* filename, MPI_Comm comm) const
+{
+#ifndef DEAL_II_COMPILER_SUPPORTS_MPI
+  //without MPI fall back to the normal way to write a vtu file:
+  (void)comm;
+
+  std::ofstream f(filename);
+  write_vtu (f);
+#else
+
+  int myrank, nproc;
+  MPI_Comm_rank(comm, &myrank);
+  MPI_Comm_size(comm, &nproc);
+
+  MPI_Info info;
+  MPI_Info_create(&info);
+  MPI_File fh;
+  MPI_File_open(MPI_COMM_WORLD, const_cast<char*>(filename),
+                MPI_MODE_CREATE | MPI_MODE_WRONLY, info, &fh);
+  MPI_File_set_size(fh, 0); // delete the file contents
+  // this barrier is necessary, because otherwise others might already
+  // write while one core is still setting the size to zero.
+  MPI_Barrier(comm);
+  MPI_Info_free(&info);
+
+  unsigned int header_size;
+
+  //write header
+  if (myrank==0)
+    {
+      std::stringstream ss;
+      DataOutBase::write_vtu_header(ss);
+      header_size = ss.str().size();
+      MPI_File_write(fh, const_cast<char*>(ss.str().c_str()), header_size, MPI_CHAR, NULL);
+    }
+
+  MPI_Bcast(&header_size, 1, MPI_INT, 0, comm);
+
+  MPI_File_seek_shared( fh, header_size, MPI_SEEK_SET );
+  {
+    std::stringstream ss;
+    DataOutBase::write_vtu_main (get_patches(), get_dataset_names(),
+                                 get_vector_data_ranges(),
+                                 vtk_flags, ss);
+    MPI_File_write_ordered(fh, const_cast<char*>(ss.str().c_str()), ss.str().size(), MPI_CHAR, NULL);
+  }
+
+  //write footer
+  if (myrank==0)
+    {
+      std::stringstream ss;
+      DataOutBase::write_vtu_footer(ss);
+      unsigned int footer_size = ss.str().size();
+      MPI_File_write_shared(fh, const_cast<char*>(ss.str().c_str()), footer_size, MPI_CHAR, NULL);
+    }
+  MPI_File_close( &fh );
+#endif
+}
+
 
 
 template <int dim, int spacedim>

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.