[cig-commits] commit by bangerth to /var/svn/dealii/aspect
dealii.demon at gmail.com
dealii.demon at gmail.com
Tue Aug 13 08:30:37 PDT 2013
Revision 1838
Move tracer output classes into their own namespace.
U trunk/aspect/include/aspect/particle/output.h
U trunk/aspect/include/aspect/particle/world.h
U trunk/aspect/include/aspect/postprocess/tracer.h
U trunk/aspect/source/particle/output.cc
U trunk/aspect/source/postprocess/tracer.cc
http://www.dealii.org/websvn/revision.php?repname=Aspect+Repository&path=%2F&rev=1838&peg=1838
Diff:
Modified: trunk/aspect/include/aspect/particle/output.h
===================================================================
--- trunk/aspect/include/aspect/particle/output.h 2013-08-13 15:22:12 UTC (rev 1837)
+++ trunk/aspect/include/aspect/particle/output.h 2013-08-13 15:30:01 UTC (rev 1838)
@@ -30,76 +30,79 @@
{
namespace Particle
{
- /**
- * Abstract base class used for classes that generate particle output
- */
- template <int dim, class T>
- class Output
+ namespace Output
{
- protected:
- /**
- * MPI communicator to be used for output synchronization
- */
- MPI_Comm communicator;
+ /**
+ * Abstract base class used for classes that generate particle output
+ */
+ template <int dim, class T>
+ class Interface
+ {
+ protected:
+ /**
+ * MPI communicator to be used for output synchronization
+ */
+ MPI_Comm communicator;
- /**
- * Internal index of file output number, must be incremented
- * by derived classes when they create a new file.
- */
- unsigned int file_index;
+ /**
+ * Internal index of file output number, must be incremented
+ * by derived classes when they create a new file.
+ */
+ unsigned int file_index;
- // Path to directory in which to put particle output files
- std::string output_dir;
+ // Path to directory in which to put particle output files
+ std::string output_dir;
- public:
- /**
- * Constructor.
- */
- Output()
- :
- file_index (0)
- {}
+ public:
+ /**
+ * Constructor.
+ */
+ Interface()
+ :
+ file_index (0)
+ {}
- unsigned int self_rank()
- {
- return Utilities::MPI::this_mpi_process(communicator);
- };
+ unsigned int self_rank()
+ {
+ return Utilities::MPI::this_mpi_process(communicator);
+ };
- unsigned int world_size()
- {
- return Utilities::MPI::n_mpi_processes(communicator);
- };
+ unsigned int world_size()
+ {
+ return Utilities::MPI::n_mpi_processes(communicator);
+ };
- void set_mpi_comm(MPI_Comm new_comm_world)
- {
- communicator = new_comm_world;
- };
+ void set_mpi_comm(MPI_Comm new_comm_world)
+ {
+ communicator = new_comm_world;
+ };
- void set_output_directory(const std::string &new_out_dir)
- {
- output_dir = new_out_dir;
- };
+ void set_output_directory(const std::string &new_out_dir)
+ {
+ output_dir = new_out_dir;
+ };
- virtual std::string output_particle_data(const std::multimap<LevelInd, T> &particles,
- const double ¤t_time) = 0;
- };
+ virtual std::string output_particle_data(const std::multimap<LevelInd, T> &particles,
+ const double ¤t_time) = 0;
+ };
- /**
- * Create an output object given the specified name.
- */
- template <int dim, class T>
- Output<dim, T> *
- create_output_object (const std::string &data_format_name);
+ /**
+ * Create an output object given the specified name.
+ */
+ template <int dim, class T>
+ Interface<dim, T> *
+ create_output_object (const std::string &data_format_name);
- /**
- * Return a list of names (separated by '|') of possible writers of graphical
- * output formats for particle data.
- */
- std::string
- output_object_names ();
+ /**
+ * Return a list of names (separated by '|') of possible writers of graphical
+ * output formats for particle data.
+ */
+ std::string
+ output_object_names ();
+ }
}
}
Modified: trunk/aspect/include/aspect/particle/world.h
===================================================================
--- trunk/aspect/include/aspect/particle/world.h 2013-08-13 15:22:12 UTC (rev 1837)
+++ trunk/aspect/include/aspect/particle/world.h 2013-08-13 15:30:01 UTC (rev 1838)
@@ -98,7 +98,7 @@
// Create the roulette wheel based on volumes of local cells
total_volume = 0;
for (typename parallel::distributed::Triangulation<dim>::active_cell_iterator
- it=triangulation->begin_active(); it!=triangulation->end(); ++it)
+ it=triangulation->begin_active(); it!=triangulation->end(); ++it)
{
if (it->is_locally_owned())
{
@@ -338,7 +338,7 @@
// Calculate the number of particles in this domain as a fraction of total volume
total_volume = local_volume = 0;
for (typename parallel::distributed::Triangulation<dim>::active_cell_iterator
- it=triangulation->begin_active();
+ it=triangulation->begin_active();
it!=triangulation->end(); ++it)
{
double cell_volume = it->measure();
Modified: trunk/aspect/include/aspect/postprocess/tracer.h
===================================================================
--- trunk/aspect/include/aspect/postprocess/tracer.h 2013-08-13 15:22:12 UTC (rev 1837)
+++ trunk/aspect/include/aspect/postprocess/tracer.h 2013-08-13 15:30:01 UTC (rev 1838)
@@ -46,9 +46,9 @@
Particle::Integrator<dim, Particle::BaseParticle<dim> > *integrator;
/**
- * Abstract output object
+ * Pointer to an output object
*/
- Particle::Output<dim, Particle::BaseParticle<dim> > *output;
+ Particle::Output::Interface<dim, Particle::BaseParticle<dim> > *output;
/**
* Whether this set has been initialized yet or not
Modified: trunk/aspect/source/particle/output.cc
===================================================================
--- trunk/aspect/source/particle/output.cc 2013-08-13 15:22:12 UTC (rev 1837)
+++ trunk/aspect/source/particle/output.cc 2013-08-13 15:30:01 UTC (rev 1838)
@@ -31,464 +31,467 @@
{
namespace Particle
{
- // Blank class to avoid output of data, for example, if particles are used to represent
- // other physical phenomenon that influences the simulation and we don't care about their positions
- template <int dim, class T>
- class NullOutput : public Output<dim, T>
+ namespace Output
{
- public:
- virtual std::string output_particle_data(const std::multimap<LevelInd, T> &particles, const double ¤t_time)
- {
- return "";
- };
- };
+ // Blank class to avoid output of data, for example, if particles are used to represent
+ // other physical phenomenon that influences the simulation and we don't care about their positions
+ template <int dim, class T>
+ class NullOutput : public Interface<dim, T>
+ {
+ public:
+ virtual std::string output_particle_data(const std::multimap<LevelInd, T> &particles, const double ¤t_time)
+ {
+ return "";
+ };
+ };
- template <int dim, class T>
- class ASCIIOutput : public Output<dim, T>
- {
- public:
- virtual std::string output_particle_data(const std::multimap<LevelInd, T> &particles,
- const double ¤t_time)
- {
- typename std::multimap<LevelInd, T>::const_iterator it;
- unsigned int i;
- std::string output_file_prefix, output_path_prefix, full_filename;
- std::vector<MPIDataInfo> data_info;
- std::vector<MPIDataInfo>::iterator dit;
- char *particle_data, *p;
+ template <int dim, class T>
+ class ASCIIOutput : public Interface<dim, T>
+ {
+ public:
+ virtual std::string output_particle_data(const std::multimap<LevelInd, T> &particles,
+ const double ¤t_time)
+ {
+ typename std::multimap<LevelInd, T>::const_iterator it;
+ unsigned int i;
+ std::string output_file_prefix, output_path_prefix, full_filename;
+ std::vector<MPIDataInfo> data_info;
+ std::vector<MPIDataInfo>::iterator dit;
+ char *particle_data, *p;
- output_file_prefix = "particle-" + Utilities::int_to_string (this->file_index, 5);
- output_path_prefix = this->output_dir + output_file_prefix;
- full_filename = output_path_prefix + "." + Utilities::int_to_string(this->self_rank(), 4) + ".txt";
- std::ofstream output (full_filename.c_str());
- if (!output)
- std::cout << "ERROR: proc " << this->self_rank() << " could not create " << full_filename << std::endl;
+ output_file_prefix = "particle-" + Utilities::int_to_string (this->file_index, 5);
+ output_path_prefix = this->output_dir + output_file_prefix;
+ full_filename = output_path_prefix + "." + Utilities::int_to_string(this->self_rank(), 4) + ".txt";
+ std::ofstream output (full_filename.c_str());
+ if (!output)
+ std::cout << "ERROR: proc " << this->self_rank() << " could not create " << full_filename << std::endl;
- // Get the data types
- T::add_mpi_types(data_info);
+ // Get the data types
+ T::add_mpi_types(data_info);
- // Print the header line
- output << "# ";
- for (dit=data_info.begin(); dit!=data_info.end(); ++dit)
- {
- // If it's a 1D element, print just the name, otherwise use []
- if (dit->_num_elems == 1)
- {
- output << dit->_name << " ";
- }
- else
- {
- for (i=0; i<dit->_num_elems; ++i) output << dit->_name << "[" << i << "] ";
- }
- }
- output << "
";
+ // Print the header line
+ output << "# ";
+ for (dit=data_info.begin(); dit!=data_info.end(); ++dit)
+ {
+ // If it's a 1D element, print just the name, otherwise use []
+ if (dit->_num_elems == 1)
+ {
+ output << dit->_name << " ";
+ }
+ else
+ {
+ for (i=0; i<dit->_num_elems; ++i) output << dit->_name << "[" << i << "] ";
+ }
+ }
+ output << "
";
- // And print the data for each particle
- particle_data = new char[T::data_len(HDF5_DATA)];
- for (it=particles.begin(); it!=particles.end(); ++it)
- {
- it->second.write_data(HDF5_DATA, particle_data);
- p = particle_data;
- for (dit=data_info.begin(); dit!=data_info.end(); ++dit)
- {
- // Currently assumes all data is double, may need to change this
- for (i=0; i<dit->_num_elems; ++i)
- {
- output << ((double *)p)[0] << " ";
- p += sizeof(double);
- }
- }
- output << "
";
- }
- delete[] particle_data;
+ // And print the data for each particle
+ particle_data = new char[T::data_len(HDF5_DATA)];
+ for (it=particles.begin(); it!=particles.end(); ++it)
+ {
+ it->second.write_data(HDF5_DATA, particle_data);
+ p = particle_data;
+ for (dit=data_info.begin(); dit!=data_info.end(); ++dit)
+ {
+ // Currently assumes all data is double, may need to change this
+ for (i=0; i<dit->_num_elems; ++i)
+ {
+ output << ((double *)p)[0] << " ";
+ p += sizeof(double);
+ }
+ }
+ output << "
";
+ }
+ delete[] particle_data;
- output.close();
+ output.close();
- this->file_index++;
+ this->file_index++;
- return output_path_prefix;
- };
- };
+ return output_path_prefix;
+ };
+ };
- template <int dim, class T>
- class VTUOutput : public Output<dim, T>
- {
- /**
- * A list of pairs (time, pvtu_filename) that have so far been written
- * and that we will pass to DataOutInterface::write_pvd_record
- * to create a master file that can make the association
- * between simulation time and corresponding file name (this
- * is done because there is no way to store the simulation
- * time inside the .pvtu or .vtu files).
- */
+ template <int dim, class T>
+ class VTUOutput : public Interface<dim, T>
+ {
+ /**
+ * A list of pairs (time, pvtu_filename) that have so far been written
+ * and that we will pass to DataOutInterface::write_pvd_record
+ * to create a master file that can make the association
+ * between simulation time and corresponding file name (this
+ * is done because there is no way to store the simulation
+ * time inside the .pvtu or .vtu files).
+ */
//TODO: This needs to be serialized
- std::vector<std::pair<double,std::string> > times_and_pvtu_file_names;
- std::vector<std::string> vtu_file_names;
+ std::vector<std::pair<double,std::string> > times_and_pvtu_file_names;
+ std::vector<std::string> vtu_file_names;
- public:
- virtual std::string output_particle_data(const std::multimap<LevelInd, T> &particles, const double ¤t_time)
- {
- typename std::multimap<LevelInd, T>::const_iterator it;
- unsigned int d, i, num_particles;
- DataOut<dim> data_out;
- std::vector<MPIDataInfo> data_info;
- std::vector<MPIDataInfo>::iterator dit;
- unsigned int data_offset;
- char *particle_data, *p;
+ public:
+ virtual std::string output_particle_data(const std::multimap<LevelInd, T> &particles, const double ¤t_time)
+ {
+ typename std::multimap<LevelInd, T>::const_iterator it;
+ unsigned int d, i, num_particles;
+ DataOut<dim> data_out;
+ std::vector<MPIDataInfo> data_info;
+ std::vector<MPIDataInfo>::iterator dit;
+ unsigned int data_offset;
+ char *particle_data, *p;
- const std::string output_file_prefix = "particles-" + Utilities::int_to_string (this->file_index, 5);
- const std::string output_path_prefix = this->output_dir + output_file_prefix;
+ const std::string output_file_prefix = "particles-" + Utilities::int_to_string (this->file_index, 5);
+ const std::string output_path_prefix = this->output_dir + output_file_prefix;
- const std::string filename = (output_file_prefix +
- "." +
- Utilities::int_to_string(this->self_rank(), 4) +
- ".vtu");
- const std::string full_filename = (this->output_dir + filename);
+ const std::string filename = (output_file_prefix +
+ "." +
+ Utilities::int_to_string(this->self_rank(), 4) +
+ ".vtu");
+ const std::string full_filename = (this->output_dir + filename);
- std::ofstream output (full_filename.c_str());
- if (!output)
- std::cout << "ERROR: proc "
- << this->self_rank() << " could not create " << filename << std::endl;
+ std::ofstream output (full_filename.c_str());
+ if (!output)
+ std::cout << "ERROR: proc "
+ << this->self_rank() << " could not create " << filename << std::endl;
- num_particles = particles.size();
+ num_particles = particles.size();
- // Write VTU file XML
- output << "<?xml version=\"1.0\"?>
";
- output << "<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">
";
- output << " <UnstructuredGrid>
";
- output << " <Piece NumberOfPoints=\"" << num_particles << "\" NumberOfCells=\"" << num_particles << "\">
";
+ // Write VTU file XML
+ output << "<?xml version=\"1.0\"?>
";
+ output << "<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">
";
+ output << " <UnstructuredGrid>
";
+ output << " <Piece NumberOfPoints=\"" << num_particles << "\" NumberOfCells=\"" << num_particles << "\">
";
- // Go through the particles on this domain and print the position of each one
- output << " <Points>
";
- output << " <DataArray name=\"Position\" type=\"Float64\" NumberOfComponents=\"3\" Format=\"ascii\">
";
- for (it=particles.begin(); it!=particles.end(); ++it)
- {
- output << " " << it->second.location();
+ // Go through the particles on this domain and print the position of each one
+ output << " <Points>
";
+ output << " <DataArray name=\"Position\" type=\"Float64\" NumberOfComponents=\"3\" Format=\"ascii\">
";
+ for (it=particles.begin(); it!=particles.end(); ++it)
+ {
+ output << " " << it->second.location();
- // pad with zeros since VTU format wants x/y/z coordinates
- for (d=dim; d<3; ++d)
- output << " 0.0";
+ // pad with zeros since VTU format wants x/y/z coordinates
+ for (d=dim; d<3; ++d)
+ output << " 0.0";
- output << "
";
- }
- output << " </DataArray>
";
- output << " </Points>
";
+ output << "
";
+ }
+ output << " </DataArray>
";
+ output << " </Points>
";
- // Write cell related data (empty)
- output << " <Cells>
";
- output << " <DataArray type=\"Int32\" Name=\"connectivity\" Format=\"ascii\">
";
- for (i=1; i<=num_particles; ++i)
- output << " " << i-1 << "
";
- output << " </DataArray>
";
- output << " <DataArray type=\"Int32\" Name=\"offsets\" Format=\"ascii\">
";
- for (i=1; i<=num_particles; ++i)
- output << " " << i << "
";
- output << " </DataArray>
";
- output << " <DataArray type=\"UInt8\" Name=\"types\" Format=\"ascii\">
";
- for (i=1; i<=num_particles; ++i)
- output << " 1
";
- output << " </DataArray>
";
- output << " </Cells>
";
+ // Write cell related data (empty)
+ output << " <Cells>
";
+ output << " <DataArray type=\"Int32\" Name=\"connectivity\" Format=\"ascii\">
";
+ for (i=1; i<=num_particles; ++i)
+ output << " " << i-1 << "
";
+ output << " </DataArray>
";
+ output << " <DataArray type=\"Int32\" Name=\"offsets\" Format=\"ascii\">
";
+ for (i=1; i<=num_particles; ++i)
+ output << " " << i << "
";
+ output << " </DataArray>
";
+ output << " <DataArray type=\"UInt8\" Name=\"types\" Format=\"ascii\">
";
+ for (i=1; i<=num_particles; ++i)
+ output << " 1
";
+ output << " </DataArray>
";
+ output << " </Cells>
";
- // Get the data types
- T::add_mpi_types(data_info);
+ // Get the data types
+ T::add_mpi_types(data_info);
- // Write data for each particle (id, velocity, etc)
- output << " <PointData Scalars=\"scalars\">
";
+ // Write data for each particle (id, velocity, etc)
+ output << " <PointData Scalars=\"scalars\">
";
- // Print the data associated with the particles, skipping the first entry (assumed to be position)
- particle_data = new char[T::data_len(HDF5_DATA)];
- dit = data_info.begin();
- data_offset = dit->_num_elems*dit->_elem_size_bytes;
- dit++;
- for (; dit!=data_info.end(); ++dit)
- {
- output << " <DataArray type=\"Float64\" Name=\"" << dit->_name << "\" NumberOfComponents=\"" << (dit->_num_elems == 2 ? 3 : dit->_num_elems) << "\" Format=\"ascii\">
";
- for (it=particles.begin(); it!=particles.end(); ++it)
- {
- it->second.write_data(HDF5_DATA, particle_data);
- p = particle_data+data_offset;
- output << " ";
- for (d=0; d<dit->_num_elems; ++d)
- {
- output << ((double *)p)[0] << " ";
- p += sizeof(double);
- }
- if (dit->_num_elems == 2)
- output << "0 ";
- output << "
";
- }
- data_offset += dit->_num_elems*dit->_elem_size_bytes;
- output << " </DataArray>
";
- }
- delete[] particle_data;
- output << " </PointData>
";
+ // Print the data associated with the particles, skipping the first entry (assumed to be position)
+ particle_data = new char[T::data_len(HDF5_DATA)];
+ dit = data_info.begin();
+ data_offset = dit->_num_elems*dit->_elem_size_bytes;
+ dit++;
+ for (; dit!=data_info.end(); ++dit)
+ {
+ output << " <DataArray type=\"Float64\" Name=\"" << dit->_name << "\" NumberOfComponents=\"" << (dit->_num_elems == 2 ? 3 : dit->_num_elems) << "\" Format=\"ascii\">
";
+ for (it=particles.begin(); it!=particles.end(); ++it)
+ {
+ it->second.write_data(HDF5_DATA, particle_data);
+ p = particle_data+data_offset;
+ output << " ";
+ for (d=0; d<dit->_num_elems; ++d)
+ {
+ output << ((double *)p)[0] << " ";
+ p += sizeof(double);
+ }
+ if (dit->_num_elems == 2)
+ output << "0 ";
+ output << "
";
+ }
+ data_offset += dit->_num_elems*dit->_elem_size_bytes;
+ output << " </DataArray>
";
+ }
+ delete[] particle_data;
+ output << " </PointData>
";
- output << " </Piece>
";
- output << " </UnstructuredGrid>
";
- output << "</VTKFile>
";
+ output << " </Piece>
";
+ output << " </UnstructuredGrid>
";
+ output << "</VTKFile>
";
- output.close();
+ output.close();
- // Write the parallel pvtu and pvd files on the root process
- if (this->self_rank() == 0)
- {
- const std::string pvtu_filename = (output_path_prefix + ".pvtu");
+ // Write the parallel pvtu and pvd files on the root process
+ if (this->self_rank() == 0)
+ {
+ const std::string pvtu_filename = (output_path_prefix + ".pvtu");
- std::ofstream pvtu_output (pvtu_filename.c_str());
- if (!pvtu_output) std::cout << "ERROR: Could not create " << filename << std::endl;
+ std::ofstream pvtu_output (pvtu_filename.c_str());
+ if (!pvtu_output) std::cout << "ERROR: Could not create " << filename << std::endl;
- pvtu_output << "<?xml version=\"1.0\"?>
";
- pvtu_output << "<VTKFile type=\"PUnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">
";
- pvtu_output << " <PUnstructuredGrid GhostLevel=\"0\">
";
- pvtu_output << " <PPoints>
";
- pvtu_output << " <PDataArray type=\"Float64\" NumberOfComponents=\"3\" format=\"ascii\"/>
";
- pvtu_output << " </PPoints>
";
- pvtu_output << " <PPointData Scalars=\"scalars\">
";
+ pvtu_output << "<?xml version=\"1.0\"?>
";
+ pvtu_output << "<VTKFile type=\"PUnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">
";
+ pvtu_output << " <PUnstructuredGrid GhostLevel=\"0\">
";
+ pvtu_output << " <PPoints>
";
+ pvtu_output << " <PDataArray type=\"Float64\" NumberOfComponents=\"3\" format=\"ascii\"/>
";
+ pvtu_output << " </PPoints>
";
+ pvtu_output << " <PPointData Scalars=\"scalars\">
";
- for (dit=data_info.begin()+1; dit!=data_info.end(); ++dit)
- {
- pvtu_output << " <PDataArray type=\"Float64\" Name=\"" << dit->_name << "\" NumberOfComponents=\"" << (dit->_num_elems == 2 ? 3 : dit->_num_elems) << "\" format=\"ascii\"/>
";
- }
- pvtu_output << " </PPointData>
";
- for (i=0; i<this->world_size(); ++i)
- {
- pvtu_output << " <Piece Source=\"" << output_file_prefix << "." << Utilities::int_to_string(i, 4) << ".vtu\"/>
";
- }
- pvtu_output << " </PUnstructuredGrid>
";
- pvtu_output << "</VTKFile>
";
- pvtu_output.close();
+ for (dit=data_info.begin()+1; dit!=data_info.end(); ++dit)
+ {
+ pvtu_output << " <PDataArray type=\"Float64\" Name=\"" << dit->_name << "\" NumberOfComponents=\"" << (dit->_num_elems == 2 ? 3 : dit->_num_elems) << "\" format=\"ascii\"/>
";
+ }
+ pvtu_output << " </PPointData>
";
+ for (i=0; i<this->world_size(); ++i)
+ {
+ pvtu_output << " <Piece Source=\"" << output_file_prefix << "." << Utilities::int_to_string(i, 4) << ".vtu\"/>
";
+ }
+ pvtu_output << " </PUnstructuredGrid>
";
+ pvtu_output << "</VTKFile>
";
+ pvtu_output.close();
- times_and_pvtu_file_names.push_back(std::make_pair(current_time,
- pvtu_filename));
- vtu_file_names.push_back (full_filename);
+ times_and_pvtu_file_names.push_back(std::make_pair(current_time,
+ pvtu_filename));
+ vtu_file_names.push_back (full_filename);
- // write .pvd and .visit records
- const std::string pvd_master_filename = (this->output_dir + "particles.pvd");
- std::ofstream pvd_master (pvd_master_filename.c_str());
- data_out.write_pvd_record (pvd_master, times_and_pvtu_file_names);
+ // write .pvd and .visit records
+ const std::string pvd_master_filename = (this->output_dir + "particles.pvd");
+ std::ofstream pvd_master (pvd_master_filename.c_str());
+ data_out.write_pvd_record (pvd_master, times_and_pvtu_file_names);
//TODO: write a global .visit record. this needs a variant of the write_visit_record
// function
-/*
- const std::string visit_master_filename = (this->output_dir + "particles.visit");
- std::ofstream visit_master (visit_master_filename.c_str());
- data_out.write_visit_record (visit_master, vtu_file_names);
-*/
- }
- this->file_index++;
+ /*
+ const std::string visit_master_filename = (this->output_dir + "particles.visit");
+ std::ofstream visit_master (visit_master_filename.c_str());
+ data_out.write_visit_record (visit_master, vtu_file_names);
+ */
+ }
+ this->file_index++;
- return output_path_prefix;
- };
- };
+ return output_path_prefix;
+ };
+ };
- template <int dim, class T>
- class HDF5Output : public Output<dim, T>
- {
- private:
- // A set of XDMF data objects to create the XDMF file for particles
- std::vector<XDMFEntry> xdmf_entries;
+ template <int dim, class T>
+ class HDF5Output : public Interface<dim, T>
+ {
+ private:
+ // A set of XDMF data objects to create the XDMF file for particles
+ std::vector<XDMFEntry> xdmf_entries;
- public:
- virtual std::string output_particle_data(const std::multimap<LevelInd, T> &particles, const double ¤t_time)
- {
- typename std::multimap<LevelInd, T>::const_iterator it;
- std::string output_file_prefix, output_path_prefix, full_filename;
+ public:
+ virtual std::string output_particle_data(const std::multimap<LevelInd, T> &particles, const double ¤t_time)
+ {
+ typename std::multimap<LevelInd, T>::const_iterator it;
+ std::string output_file_prefix, output_path_prefix, full_filename;
- output_file_prefix = "particle-" + Utilities::int_to_string (this->file_index, 5);
- output_path_prefix = this->output_dir + output_file_prefix;
+ output_file_prefix = "particle-" + Utilities::int_to_string (this->file_index, 5);
+ output_path_prefix = this->output_dir + output_file_prefix;
#ifdef DEAL_II_HAVE_HDF5
- // TODO: error checking for H5 calls
- hid_t h5_file_id;
- hid_t plist_id, xfer_plist_id, position_dataset_id, velocity_dataset_id, pid_dataset_id;
- hid_t file_dataspace_id, pos_file_dataspace_id, vel_file_dataspace_id, pid_file_dataspace_id;
- hid_t dim_dataspace_id, one_dim_ds_id;
- hid_t dim_mem_ds_id, one_dim_mem_ds_id;
- hid_t pattr_id;
- herr_t status;
- hsize_t dims[2], offset[2], count[2];
- unsigned int mpi_offset, mpi_count, local_particle_count, global_particle_count, d, i;
+ // TODO: error checking for H5 calls
+ hid_t h5_file_id;
+ hid_t plist_id, xfer_plist_id, position_dataset_id, velocity_dataset_id, pid_dataset_id;
+ hid_t file_dataspace_id, pos_file_dataspace_id, vel_file_dataspace_id, pid_file_dataspace_id;
+ hid_t dim_dataspace_id, one_dim_ds_id;
+ hid_t dim_mem_ds_id, one_dim_mem_ds_id;
+ hid_t pattr_id;
+ herr_t status;
+ hsize_t dims[2], offset[2], count[2];
+ unsigned int mpi_offset, mpi_count, local_particle_count, global_particle_count, d, i;
- double *pos_data, *vel_data, *id_data;
+ double *pos_data, *vel_data, *id_data;
- std::string h5_filename = output_path_prefix+".h5";
+ std::string h5_filename = output_path_prefix+".h5";
- // Create parallel file access
- plist_id = H5Pcreate(H5P_FILE_ACCESS);
+ // Create parallel file access
+ plist_id = H5Pcreate(H5P_FILE_ACCESS);
#ifdef H5_HAVE_PARALLEL
- H5Pset_fapl_mpio(plist_id, this->communicator, MPI_INFO_NULL);
+ H5Pset_fapl_mpio(plist_id, this->communicator, MPI_INFO_NULL);
#endif
- // Create the file
- h5_file_id = H5Fcreate(h5_filename.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, plist_id);
+ // Create the file
+ h5_file_id = H5Fcreate(h5_filename.c_str(), H5F_ACC_TRUNC, H5P_DEFAULT, plist_id);
- // Create the file dataspace descriptions
- local_particle_count = particles.size();
- // TODO: error checking on MPI call
- MPI_Comm com = this->communicator;
- MPI_Allreduce(&local_particle_count, &global_particle_count, 1, MPI_UNSIGNED, MPI_SUM, com);
- dims[0] = global_particle_count;
- dims[1] = 3;
- one_dim_ds_id = H5Screate_simple(1, dims, NULL);
- dim_dataspace_id = H5Screate_simple(2, dims, NULL);
+ // Create the file dataspace descriptions
+ local_particle_count = particles.size();
+ // TODO: error checking on MPI call
+ MPI_Comm com = this->communicator;
+ MPI_Allreduce(&local_particle_count, &global_particle_count, 1, MPI_UNSIGNED, MPI_SUM, com);
+ dims[0] = global_particle_count;
+ dims[1] = 3;
+ one_dim_ds_id = H5Screate_simple(1, dims, NULL);
+ dim_dataspace_id = H5Screate_simple(2, dims, NULL);
- // Create the datasets
+ // Create the datasets
#if H5Dcreate_vers == 1
- position_dataset_id = H5Dcreate(h5_file_id, "nodes", H5T_NATIVE_DOUBLE, dim_dataspace_id, H5P_DEFAULT);
- velocity_dataset_id = H5Dcreate(h5_file_id, "velocity", H5T_NATIVE_DOUBLE, dim_dataspace_id, H5P_DEFAULT);
- pid_dataset_id = H5Dcreate(h5_file_id, "id", H5T_NATIVE_DOUBLE, one_dim_ds_id, H5P_DEFAULT);
+ position_dataset_id = H5Dcreate(h5_file_id, "nodes", H5T_NATIVE_DOUBLE, dim_dataspace_id, H5P_DEFAULT);
+ velocity_dataset_id = H5Dcreate(h5_file_id, "velocity", H5T_NATIVE_DOUBLE, dim_dataspace_id, H5P_DEFAULT);
+ pid_dataset_id = H5Dcreate(h5_file_id, "id", H5T_NATIVE_DOUBLE, one_dim_ds_id, H5P_DEFAULT);
#else
- position_dataset_id = H5Dcreate(h5_file_id, "nodes", H5T_NATIVE_DOUBLE, dim_dataspace_id, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
- velocity_dataset_id = H5Dcreate(h5_file_id, "velocity", H5T_NATIVE_DOUBLE, dim_dataspace_id, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
- pid_dataset_id = H5Dcreate(h5_file_id, "id", H5T_NATIVE_DOUBLE, one_dim_ds_id, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
+ position_dataset_id = H5Dcreate(h5_file_id, "nodes", H5T_NATIVE_DOUBLE, dim_dataspace_id, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
+ velocity_dataset_id = H5Dcreate(h5_file_id, "velocity", H5T_NATIVE_DOUBLE, dim_dataspace_id, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
+ pid_dataset_id = H5Dcreate(h5_file_id, "id", H5T_NATIVE_DOUBLE, one_dim_ds_id, H5P_DEFAULT, H5P_DEFAULT, H5P_DEFAULT);
#endif
- // Close the file dataspaces
- H5Sclose(dim_dataspace_id);
- H5Sclose(one_dim_ds_id);
+ // Close the file dataspaces
+ H5Sclose(dim_dataspace_id);
+ H5Sclose(one_dim_ds_id);
- // Just in case we forget what they are
- count[0] = 1;
- file_dataspace_id = H5Screate_simple (1, count, NULL);
+ // Just in case we forget what they are
+ count[0] = 1;
+ file_dataspace_id = H5Screate_simple (1, count, NULL);
#if H5Acreate_vers == 1
- pattr_id = H5Acreate(h5_file_id, "Ermahgerd! Pertecrs!", H5T_NATIVE_INT, file_dataspace_id, H5P_DEFAULT);
+ pattr_id = H5Acreate(h5_file_id, "Ermahgerd! Pertecrs!", H5T_NATIVE_INT, file_dataspace_id, H5P_DEFAULT);
#else
- pattr_id = H5Acreate(h5_file_id, "Ermahgerd! Pertecrs!", H5T_NATIVE_INT, file_dataspace_id, H5P_DEFAULT, H5P_DEFAULT);
+ pattr_id = H5Acreate(h5_file_id, "Ermahgerd! Pertecrs!", H5T_NATIVE_INT, file_dataspace_id, H5P_DEFAULT, H5P_DEFAULT);
#endif
- H5Aclose(pattr_id);
- H5Sclose(file_dataspace_id);
+ H5Aclose(pattr_id);
+ H5Sclose(file_dataspace_id);
- // Count the number of particles on this process and get the offset among all processes
- mpi_count = particles.size();
- MPI_Scan(&mpi_count, &mpi_offset, 1, MPI_UNSIGNED, MPI_SUM, this->communicator);
- count[0] = mpi_count;
- count[1] = 3;
- offset[0] = mpi_offset-mpi_count;
- offset[1] = 0;
+ // Count the number of particles on this process and get the offset among all processes
+ mpi_count = particles.size();
+ MPI_Scan(&mpi_count, &mpi_offset, 1, MPI_UNSIGNED, MPI_SUM, this->communicator);
+ count[0] = mpi_count;
+ count[1] = 3;
+ offset[0] = mpi_offset-mpi_count;
+ offset[1] = 0;
- // Select the appropriate dataspace for this process
- one_dim_mem_ds_id = H5Screate_simple(1, count, NULL);
- dim_mem_ds_id = H5Screate_simple(2, count, NULL);
- pos_file_dataspace_id = H5Dget_space(position_dataset_id);
- vel_file_dataspace_id = H5Dget_space(velocity_dataset_id);
- pid_file_dataspace_id = H5Dget_space(pid_dataset_id);
+ // Select the appropriate dataspace for this process
+ one_dim_mem_ds_id = H5Screate_simple(1, count, NULL);
+ dim_mem_ds_id = H5Screate_simple(2, count, NULL);
+ pos_file_dataspace_id = H5Dget_space(position_dataset_id);
+ vel_file_dataspace_id = H5Dget_space(velocity_dataset_id);
+ pid_file_dataspace_id = H5Dget_space(pid_dataset_id);
- // And select the hyperslabs from each dataspace
- status = H5Sselect_hyperslab(pos_file_dataspace_id, H5S_SELECT_SET, offset, NULL, count, NULL);
- status = H5Sselect_hyperslab(vel_file_dataspace_id, H5S_SELECT_SET, offset, NULL, count, NULL);
- status = H5Sselect_hyperslab(pid_file_dataspace_id, H5S_SELECT_SET, offset, NULL, count, NULL);
+ // And select the hyperslabs from each dataspace
+ status = H5Sselect_hyperslab(pos_file_dataspace_id, H5S_SELECT_SET, offset, NULL, count, NULL);
+ status = H5Sselect_hyperslab(vel_file_dataspace_id, H5S_SELECT_SET, offset, NULL, count, NULL);
+ status = H5Sselect_hyperslab(pid_file_dataspace_id, H5S_SELECT_SET, offset, NULL, count, NULL);
- // Create property list for collective dataset write
- xfer_plist_id = H5Pcreate(H5P_DATASET_XFER);
+ // Create property list for collective dataset write
+ xfer_plist_id = H5Pcreate(H5P_DATASET_XFER);
#ifdef H5_HAVE_PARALLEL
- H5Pset_dxpl_mpio(xfer_plist_id, H5FD_MPIO_COLLECTIVE);
+ H5Pset_dxpl_mpio(xfer_plist_id, H5FD_MPIO_COLLECTIVE);
#endif
- // Read the local particle data
- pos_data = new double[3*particles.size()];
- vel_data = new double[3*particles.size()];
- id_data = new double[particles.size()];
+ // Read the local particle data
+ pos_data = new double[3*particles.size()];
+ vel_data = new double[3*particles.size()];
+ id_data = new double[particles.size()];
- for (i=0,it=particles.begin(); it!=particles.end(); ++i,++it)
- {
- for (d=0; d<dim; ++d)
- {
- pos_data[i*3+d] = it->second.location()(d);
- vel_data[i*3+d] = it->second.velocity()(d);
- }
- if (dim < 3)
- {
- pos_data[i*3+2] = 0;
- vel_data[i*3+2] = 0;
- }
- id_data[i] = it->second.id_num();
- }
+ for (i=0,it=particles.begin(); it!=particles.end(); ++i,++it)
+ {
+ for (d=0; d<dim; ++d)
+ {
+ pos_data[i*3+d] = it->second.location()(d);
+ vel_data[i*3+d] = it->second.velocity()(d);
+ }
+ if (dim < 3)
+ {
+ pos_data[i*3+2] = 0;
+ vel_data[i*3+2] = 0;
+ }
+ id_data[i] = it->second.id_num();
+ }
- // Write particle data to the HDF5 file
- H5Dwrite(position_dataset_id, H5T_NATIVE_DOUBLE, dim_mem_ds_id, pos_file_dataspace_id, xfer_plist_id, pos_data);
- H5Dwrite(velocity_dataset_id, H5T_NATIVE_DOUBLE, dim_mem_ds_id, vel_file_dataspace_id, xfer_plist_id, vel_data);
- H5Dwrite(pid_dataset_id, H5T_NATIVE_DOUBLE, one_dim_mem_ds_id, pid_file_dataspace_id, xfer_plist_id, id_data);
+ // Write particle data to the HDF5 file
+ H5Dwrite(position_dataset_id, H5T_NATIVE_DOUBLE, dim_mem_ds_id, pos_file_dataspace_id, xfer_plist_id, pos_data);
+ H5Dwrite(velocity_dataset_id, H5T_NATIVE_DOUBLE, dim_mem_ds_id, vel_file_dataspace_id, xfer_plist_id, vel_data);
+ H5Dwrite(pid_dataset_id, H5T_NATIVE_DOUBLE, one_dim_mem_ds_id, pid_file_dataspace_id, xfer_plist_id, id_data);
- // Clear allocated resources
- delete[] pos_data;
- delete[] vel_data;
- delete[] id_data;
- status = H5Pclose(xfer_plist_id);
- status = H5Sclose(one_dim_mem_ds_id);
- status = H5Sclose(dim_mem_ds_id);
- status = H5Sclose(pos_file_dataspace_id);
- status = H5Sclose(vel_file_dataspace_id);
- status = H5Sclose(pid_file_dataspace_id);
- status = H5Dclose(position_dataset_id);
- status = H5Dclose(velocity_dataset_id);
- status = H5Dclose(pid_dataset_id);
- status = H5Pclose(plist_id);
- status = H5Fclose(h5_file_id);
+ // Clear allocated resources
+ delete[] pos_data;
+ delete[] vel_data;
+ delete[] id_data;
+ status = H5Pclose(xfer_plist_id);
+ status = H5Sclose(one_dim_mem_ds_id);
+ status = H5Sclose(dim_mem_ds_id);
+ status = H5Sclose(pos_file_dataspace_id);
+ status = H5Sclose(vel_file_dataspace_id);
+ status = H5Sclose(pid_file_dataspace_id);
+ status = H5Dclose(position_dataset_id);
+ status = H5Dclose(velocity_dataset_id);
+ status = H5Dclose(pid_dataset_id);
+ status = H5Pclose(plist_id);
+ status = H5Fclose(h5_file_id);
- // Record and output XDMF info on root process
- if (this->self_rank() == 0)
- {
- std::string local_h5_filename = output_file_prefix+".h5";
- XDMFEntry entry(local_h5_filename, current_time, global_particle_count, 0, 3);
- DataOut<dim> data_out;
- const std::string xdmf_filename = (this->output_dir + "particle.xdmf");
+ // Record and output XDMF info on root process
+ if (this->self_rank() == 0)
+ {
+ std::string local_h5_filename = output_file_prefix+".h5";
+ XDMFEntry entry(local_h5_filename, current_time, global_particle_count, 0, 3);
+ DataOut<dim> data_out;
+ const std::string xdmf_filename = (this->output_dir + "particle.xdmf");
- entry.add_attribute("velocity", 3);
- entry.add_attribute("id", 1);
- xdmf_entries.push_back(entry);
+ entry.add_attribute("velocity", 3);
+ entry.add_attribute("id", 1);
+ xdmf_entries.push_back(entry);
- data_out.write_xdmf_file(xdmf_entries, xdmf_filename.c_str(), this->communicator);
- }
+ data_out.write_xdmf_file(xdmf_entries, xdmf_filename.c_str(), this->communicator);
+ }
#endif
- this->file_index++;
+ this->file_index++;
- return output_path_prefix;
- };
- };
+ return output_path_prefix;
+ };
+ };
- template <int dim, class T>
- Output<dim,T> *
- create_output_object (const std::string &data_output_format)
- {
- if (data_output_format == "ascii")
- return new Particle::ASCIIOutput<dim,T>();
- else if (data_output_format == "vtu")
- return new Particle::VTUOutput<dim,T>();
- else if (data_output_format == "hdf5")
- return new Particle::HDF5Output<dim,T>();
- else if (data_output_format == "none")
- return new Particle::NullOutput<dim,T>();
- else
- Assert (false, ExcNotImplemented());
+ template <int dim, class T>
+ Interface<dim,T> *
+ create_output_object (const std::string &data_output_format)
+ {
+ if (data_output_format == "ascii")
+ return new ASCIIOutput<dim,T>();
+ else if (data_output_format == "vtu")
+ return new VTUOutput<dim,T>();
+ else if (data_output_format == "hdf5")
+ return new HDF5Output<dim,T>();
+ else if (data_output_format == "none")
+ return new NullOutput<dim,T>();
+ else
+ Assert (false, ExcNotImplemented());
- return 0;
- }
+ return 0;
+ }
- std::string
- output_object_names ()
- {
- return ("none|"
- "ascii|"
- "vtu|"
- "hdf5");
- }
More information about the CIG-COMMITS
mailing list