#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
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
write_gmv_reorder_data_vectors (const std::vector<Patch<dim,spacedim> > &patches,
Table<2,double> &data_vectors);
+
+
};
*/
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
#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>
#include <ctime>
#include <cmath>
#include <set>
-
#include <sstream>
+#include <fstream>
#ifdef HAVE_LIBZ
# include <zlib.h>
}
+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());
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";
// 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 ();
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>