From: bangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Date: Tue, 1 Jun 2010 18:44:58 +0000 (+0000)
Subject: write_pvtu_record
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=c3ea89b47e2f5a0dd8c9563932842c3d1a3b0535;p=dealii-svn.git

write_pvtu_record

git-svn-id: https://svn.dealii.org/trunk@21141 0785d39b-7218-0410-832d-ea1e28bc413d
---

diff --git a/deal.II/base/include/base/data_out_base.h b/deal.II/base/include/base/data_out_base.h
index bc2c8748ea..3310be5c76 100644
--- a/deal.II/base/include/base/data_out_base.h
+++ b/deal.II/base/include/base/data_out_base.h
@@ -1707,6 +1707,11 @@ class DataOutBase
  * provisions that allow these components to be output by a single
  * name rather than having to group several scalar fields into a
  * vector later on in the visualization program.
+ *
+ * Some visualization programs, such as ParaView, can read several separate
+ * VTU files to parallelize visualization. In that case, you need a
+ * <code>.pvtu</code> file that describes which VTU files form a group. The
+ * DataOutInterface::write_pvtu() function can generate such a master record.
  */
     template <int dim, int spacedim>
     static void write_vtu (const std::vector<Patch<dim,spacedim> > &patches,
@@ -2225,8 +2230,55 @@ class DataOutInterface : private DataOutBase
 				      * and write it to <tt>out</tt>
 				      * in Vtu (VTK's XML) format. See
 				      * DataOutBase::write_vtu.
+				      *
+				      * Some visualization programs, such as
+				      * ParaView, can read several separate
+				      * VTU files to parallelize
+				      * visualization. In that case, you need
+				      * a <code>.pvtu</code> file that
+				      * describes which VTU files form a
+				      * group. The
+				      * DataOutInterface::write_pvtu()
+				      * function can generate such a master
+				      * record.
 				      */
     void write_vtu (std::ostream &out) const;
+    
+				     /**
+				      * Some visualization programs, such as
+				      * ParaView, can read several separate
+				      * VTU files to parallelize
+				      * visualization. In that case, you need
+				      * a <code>.pvtu</code> file that
+				      * describes which VTU files (written,
+				      * for example, through the write_vtu()
+				      * function) form a group. The current
+				      * function can generate such a master
+				      * record.
+				      *
+				      * The file so written contains a list of
+				      * (scalar or vector) fields whose values
+				      * are described by the individual files
+				      * that comprise the set of parallel VTU
+				      * files along with the names of these
+				      * files. This function gets the names
+				      * and types of fields through the
+				      * get_patches() function of this class
+				      * like all the other write_xxx()
+				      * functions. The second argument to this
+				      * function specifies the names of the
+				      * files that form the parallel set.
+				      *
+				      * @note See DataOutBase::write_vtu for
+				      * writing each piece. Also note that
+				      * only one parallel process needs to
+				      * call the current function, listing the
+				      * names of the files written by all
+				      * parallel processes.
+				      */
+    void write_pvtu_record (std::ostream &out, 
+			    const std::vector<std::string> &piece_names) const;    
+    
 
     				     /**
 				      * Obtain data through get_patches()
diff --git a/deal.II/base/source/data_out_base.cc b/deal.II/base/source/data_out_base.cc
index 50f6ada17b..9f54c0ea2a 100644
--- a/deal.II/base/source/data_out_base.cc
+++ b/deal.II/base/source/data_out_base.cc
@@ -4287,15 +4287,15 @@ DataOutBase::write_vtu (const std::vector<Patch<dim,spacedim> > &patches,
       out << "    <DataArray type=\"Float32\" Name=\"";
       
       if (std_cxx1x::get<2>(vector_data_ranges[n_th_vector]) != "")
-	out << std_cxx1x::get<2>(vector_data_ranges[n_th_vector]);
+			out << std_cxx1x::get<2>(vector_data_ranges[n_th_vector]);
       else
-	{
-	  for (unsigned int i=std_cxx1x::get<0>(vector_data_ranges[n_th_vector]);
-	       i<std_cxx1x::get<1>(vector_data_ranges[n_th_vector]);
-	       ++i)
-	    out << data_names[i] << "__";
-	  out << data_names[std_cxx1x::get<1>(vector_data_ranges[n_th_vector])];
-	}
+		{
+		  for (unsigned int i=std_cxx1x::get<0>(vector_data_ranges[n_th_vector]);
+			   i<std_cxx1x::get<1>(vector_data_ranges[n_th_vector]);
+			   ++i)
+			out << data_names[i] << "__";
+		  out << data_names[std_cxx1x::get<1>(vector_data_ranges[n_th_vector])];
+		}
 
       out << "\" NumberOfComponents=\"3\" format=\"ascii\"> \n";
 
@@ -4593,6 +4593,109 @@ void DataOutInterface<dim,spacedim>::write_vtu (std::ostream &out) const
 }
 
 
+template <int dim, int spacedim>
+void
+DataOutInterface<dim,spacedim>::write_pvtu_record (std::ostream &out,
+						   const std::vector<std::string> &piece_names) const
+{
+  AssertThrow (out, ExcIO());
+	  
+  const std::vector<std::string> data_names = get_dataset_names();
+  const std::vector<std_cxx1x::tuple<unsigned int, unsigned int, std::string> > vector_data_ranges
+    = get_vector_data_ranges();
+	
+  const unsigned int n_data_sets = data_names.size();
+	
+  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=\"PUnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\"> \n";
+  out << "  <PUnstructuredGrid GhostLevel=\"0\"> \n";
+  out << "    <PPointData Scalars=\"scalars\"> \n";
+	
+				   // We need to output in the same order as the write_vtu function does:
+	
+  std::vector<bool> data_set_written (n_data_sets, false);
+  for (unsigned int n_th_vector=0; n_th_vector<vector_data_ranges.size(); ++n_th_vector)
+    {
+      AssertThrow (std_cxx1x::get<1>(vector_data_ranges[n_th_vector]) >=
+		   std_cxx1x::get<0>(vector_data_ranges[n_th_vector]),
+		   ExcLowerRange (std_cxx1x::get<1>(vector_data_ranges[n_th_vector]),
+				  std_cxx1x::get<0>(vector_data_ranges[n_th_vector])));
+      AssertThrow (std_cxx1x::get<1>(vector_data_ranges[n_th_vector]) < n_data_sets,
+		   ExcIndexRange (std_cxx1x::get<1>(vector_data_ranges[n_th_vector]),
+				  0, n_data_sets));
+      AssertThrow (std_cxx1x::get<1>(vector_data_ranges[n_th_vector]) + 1
+		   - std_cxx1x::get<0>(vector_data_ranges[n_th_vector]) <= 3,
+		   ExcMessage ("Can't declare a vector with more than 3 components "
+			       "in VTK"));
+		
+				       // mark these components as already
+				       // written:
+      for (unsigned int i=std_cxx1x::get<0>(vector_data_ranges[n_th_vector]);
+	   i<=std_cxx1x::get<1>(vector_data_ranges[n_th_vector]);
+	   ++i)
+	data_set_written[i] = true;
+		
+				       // write the
+				       // header. concatenate all the
+				       // component names with double
+				       // underscores unless a vector
+				       // name has been specified
+      out << "    <PDataArray type=\"Float32\" Name=\"";
+		
+      if (std_cxx1x::get<2>(vector_data_ranges[n_th_vector]) != "")
+	out << std_cxx1x::get<2>(vector_data_ranges[n_th_vector]);
+      else
+	{
+	  for (unsigned int i=std_cxx1x::get<0>(vector_data_ranges[n_th_vector]);
+	       i<std_cxx1x::get<1>(vector_data_ranges[n_th_vector]);
+	       ++i)
+	    out << data_names[i] << "__";
+	  out << data_names[std_cxx1x::get<1>(vector_data_ranges[n_th_vector])];
+	}
+	
+      out << "\" NumberOfComponents=\"3\" format=\"ascii\"/> \n";
+    }
+	
+  for (unsigned int data_set=0; data_set<n_data_sets; ++data_set)
+    if (data_set_written[data_set] == false)
+      {
+	out << "    <PDataArray type=\"Float32\" Name=\""
+	    << data_names[data_set]
+	    << "\" format=\"ascii\"/> \n";
+      }  
+
+  out << "    </PPointData> \n";
+
+  out << "    <PPoints> \n";
+  out << "      <PDataArray type=\"Float32\" NumberOfComponents=\"3\"/> \n";
+  out << "    </PPoints> \n";
+	
+  for(unsigned int i=0; i<piece_names.size(); ++i)		
+    out << "    <Piece Source=\"" << piece_names[i] << "\"/> \n";
+		
+  out << "  </PUnstructuredGrid> \n";
+  out << "</VTKFile> \n";
+	
+  out.flush();
+	
+				   // assert the stream is still ok
+  AssertThrow (out, ExcIO());
+}
+
+
 template <int dim, int spacedim>
 void DataOutInterface<dim,spacedim>::
 write_deal_II_intermediate (std::ostream &out) const
diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h
index e216d88a25..9d0d4e17b2 100644
--- a/deal.II/doc/news/changes.h
+++ b/deal.II/doc/news/changes.h
@@ -279,9 +279,11 @@ inconvenience this causes.
 <ol>
   <li><p>New: The DataOutBase class (and all derived classes such as DataOut,
   MatrixOut, etc) can now produce the XML-based version of the VTK file format
-  (the so-called VTU format).
+  (the so-called VTU format). Furthermore, the
+  DataOutInterfaces::write_pvtu_record function can be used to describe a set
+  of parallel VTU files as part of a single visualization set.
   <br>
-  (Scott Miller 2010/05/12)
+  (Scott Miller 2010/06/01)
   </p></li>
 
   <li><p>Changed: The Function::vector_gradient_list function was previously