From: Rene Gassmoeller Date: Tue, 17 Oct 2017 21:27:40 +0000 (-0600) Subject: Restructure DataOutFilter X-Git-Tag: v9.0.0-rc1~865^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=594ba129442bdcb4e36e63c57007ff734cd88f4f;p=dealii.git Restructure DataOutFilter --- diff --git a/include/deal.II/base/data_out_base.h b/include/deal.II/base/data_out_base.h index 125846ddbe..857c9e1c6e 100644 --- a/include/deal.II/base/data_out_base.h +++ b/include/deal.II/base/data_out_base.h @@ -1274,76 +1274,17 @@ namespace DataOutBase */ class DataOutFilter { - private: + public: /** - * Empty class to provide comparison function for Map3DPoint. + * Default constructor. */ - struct Point3Comp - { - bool operator() (const Point<3> &lhs, const Point<3> &rhs) const - { - return (lhs(0) < rhs(0) || (!(rhs(0) < lhs(0)) && (lhs(1) < rhs(1) || (!(rhs(1) < lhs(1)) && lhs(2) < rhs(2))))); - } - }; - - typedef std::multimap, unsigned int, Point3Comp> Map3DPoint; - - /// Flags used to specify filtering behavior - DataOutBase::DataOutFilterFlags flags; + DataOutFilter(); /** - * The number of space dimensions in which the vertices represented - * by the current object live. This corresponds to the usual - * @p dim argument, but since this class is not templated on the - * dimension, we need to store it here. + * Destructor with a given set of flags. See DataOutFilterFlags for + * possible flags. */ - unsigned int node_dim; - - /** - * The number of vertices per cell. Equal to - * GeometryInfo::vertices_per_cell. We need to store - * it as a run-time variable here because the dimension - * node_dim is also a run-time variable. - */ - unsigned int vertices_per_cell; - - /// Map of points to an internal index - Map3DPoint existing_points; - - /// Map of actual point index to internal point index - std::map filtered_points; - - /// Map of cells to the filtered points - std::map filtered_cells; - - /// Data set names - std::vector data_set_names; - - /// Data set dimensions - std::vector data_set_dims; - - /// Data set data - std::vector > data_sets; - - /** - * Record a cell vertex index based on the internal reordering. - */ - void internal_add_cell(const unsigned int &cell_index, const unsigned int &pt_index); - - public: - DataOutFilter() - : - flags(false, true), - node_dim (numbers::invalid_unsigned_int), - vertices_per_cell (numbers::invalid_unsigned_int) - {} - - DataOutFilter(const DataOutBase::DataOutFilterFlags &flags) - : - flags(flags), - node_dim (numbers::invalid_unsigned_int), - vertices_per_cell (numbers::invalid_unsigned_int) - {} + DataOutFilter(const DataOutBase::DataOutFilterFlags &flags); /** * Write a point with the specified index into the filtered data set. If @@ -1351,13 +1292,18 @@ namespace DataOutBase * provided index will internally refer to another recorded point. */ template - void write_point(const unsigned int &index, const Point &p); + void write_point(const unsigned int index, + const Point &p); /** * Record a deal.II cell in the internal reordered format. */ template - void write_cell(unsigned int index, unsigned int start, unsigned int d1, unsigned int d2, unsigned int d3); + void write_cell(const unsigned int index, + const unsigned int start, + const unsigned int d1, + const unsigned int d2, + const unsigned int d3); /** * Filter and record a data set. If there are multiple values at a given @@ -1365,7 +1311,10 @@ namespace DataOutBase * chosen as the recorded value. In the future this can be expanded to * average/min/max multiple values at a given vertex. */ - void write_data_set(const std::string &name, const unsigned int &dimension, const unsigned int &set_num, const Table<2,double> &data_vectors); + void write_data_set(const std::string &name, + const unsigned int dimension, + const unsigned int set_num, + const Table<2,double> &data_vectors); /** * Resize and fill a vector with all the filtered node vertex points, for @@ -1377,69 +1326,138 @@ namespace DataOutBase * Resize and fill a vector with all the filtered cell vertex indices, for * output to a file. */ - void fill_cell_data(const unsigned int &local_node_offset, std::vector &cell_data) const; + void fill_cell_data(const unsigned int local_node_offset, + std::vector &cell_data) const; /** * Get the name of the data set indicated by the set number. */ - std::string get_data_set_name(const unsigned int &set_num) const - { - return data_set_names.at(set_num); - }; + std::string get_data_set_name(const unsigned int set_num) const; /** * Get the dimensionality of the data set indicated by the set number. */ - unsigned int get_data_set_dim(const unsigned int &set_num) const - { - return data_set_dims.at(set_num); - }; + unsigned int get_data_set_dim(const unsigned int set_num) const; /** * Get the raw double valued data of the data set indicated by the set * number. */ - const double *get_data_set(const unsigned int &set_num) const - { - return &data_sets[set_num][0]; - }; + const double *get_data_set(const unsigned int set_num) const; /** * Return the number of nodes in this DataOutFilter. This may be smaller * than the original number of nodes if filtering is enabled. */ - unsigned int n_nodes() const - { - return existing_points.size(); - }; + unsigned int n_nodes() const; /** * Return the number of filtered cells in this DataOutFilter. Cells are * not filtered so this will be the original number of cells. */ - unsigned int n_cells() const - { - return filtered_cells.size()/vertices_per_cell; - }; + unsigned int n_cells() const; /** * Return the number of filtered data sets in this DataOutFilter. Data * sets are not filtered so this will be the original number of data sets. */ - unsigned int n_data_sets() const - { - return data_set_names.size(); - }; + unsigned int n_data_sets() const; /** * Empty functions to do base class inheritance. */ - void flush_points () {}; + void flush_points (); /** * Empty functions to do base class inheritance. */ - void flush_cells () {}; + void flush_cells (); + + + private: + /** + * Empty class to provide comparison function for Map3DPoint. + */ + struct Point3Comp + { + bool operator() (const Point<3> &one, + const Point<3> &two) const + { + /* + * The return statement below is an optimized version of the following code: + * + * for (unsigned int d=0; d<3; ++d) + * { + * if (one(d) < two(d)) + * return true; + * else if (one(d) > two(d)) + * return false; + * } + * return false; + */ + + return (one(0) < two(0) || (!(two(0) < one(0)) && (one(1) < two(1) || (!(two(1) < one(1)) && one(2) < two(2))))); + } + }; + + typedef std::multimap, unsigned int, Point3Comp> Map3DPoint; + + /** + * Flags used to specify filtering behavior. + */ + DataOutBase::DataOutFilterFlags flags; + + /** + * The number of space dimensions in which the vertices represented + * by the current object live. This corresponds to the usual + * @p dim argument, but since this class is not templated on the + * dimension, we need to store it here. + */ + unsigned int node_dim; + + /** + * The number of vertices per cell. Equal to + * GeometryInfo::vertices_per_cell. We need to store + * it as a run-time variable here because the dimension + * node_dim is also a run-time variable. + */ + unsigned int vertices_per_cell; + + /** + * Map of points to an internal index. + */ + Map3DPoint existing_points; + + /** + * Map of actual point index to internal point index. + */ + std::map filtered_points; + + /** + * Map of cells to the filtered points. + */ + std::map filtered_cells; + + /** + * Data set names. + */ + std::vector data_set_names; + + /** + * Data set dimensions. + */ + std::vector data_set_dims; + + /** + * Data set data. + */ + std::vector > data_sets; + + /** + * Record a cell vertex index based on the internal reordering. + */ + void internal_add_cell(const unsigned int cell_index, + const unsigned int pt_index); }; diff --git a/source/base/data_out_base.cc b/source/base/data_out_base.cc index 74af2846b9..4200ca4441 100644 --- a/source/base/data_out_base.cc +++ b/source/base/data_out_base.cc @@ -453,124 +453,221 @@ namespace DataOutBase ExcInternalError()); } } -} -//----------------------------------------------------------------------// -// DataOutFilter class member functions -//----------------------------------------------------------------------// -template -void DataOutBase::DataOutFilter::write_point(const unsigned int &index, const Point &p) -{ - Map3DPoint::const_iterator it; - unsigned int internal_ind; - Point<3> int_pt; - for (int d=0; d<3; ++d) int_pt(d) = (d < dim ? p(d) : 0); - node_dim = dim; - it = existing_points.find(int_pt); + DataOutFilter::DataOutFilter() + : + flags(false, true), + node_dim (numbers::invalid_unsigned_int), + vertices_per_cell (numbers::invalid_unsigned_int) + {} - // If the point isn't in the set, or we're not filtering duplicate points, add it - if (it == existing_points.end() || !flags.filter_duplicate_vertices) - { - internal_ind = existing_points.size(); - existing_points.insert(std::make_pair(int_pt, internal_ind)); - } - else - { - internal_ind = it->second; - } - // Now add the index to the list of filtered points - filtered_points[index] = internal_ind; -} -void DataOutBase::DataOutFilter::internal_add_cell(const unsigned int &cell_index, const unsigned int &pt_index) -{ - filtered_cells[cell_index] = filtered_points[pt_index]; -} -void DataOutBase::DataOutFilter::fill_node_data(std::vector &node_data) const -{ - Map3DPoint::const_iterator it; + DataOutFilter::DataOutFilter(const DataOutBase::DataOutFilterFlags &flags) + : + flags(flags), + node_dim (numbers::invalid_unsigned_int), + vertices_per_cell (numbers::invalid_unsigned_int) + {} - node_data.resize(existing_points.size()*node_dim); - for (it=existing_points.begin(); it!=existing_points.end(); ++it) - { - for (unsigned int d=0; dsecond+d] = it->first(d); - } -} -void DataOutBase::DataOutFilter::fill_cell_data(const unsigned int &local_node_offset, std::vector &cell_data) const -{ - std::map::const_iterator it; + template + void + DataOutFilter::write_point(const unsigned int index, + const Point &p) + { + node_dim = dim; - cell_data.resize(filtered_cells.size()); + Point<3> int_pt; + for (unsigned int d=0; dfirst] = it->second+local_node_offset; - } -} + const Map3DPoint::const_iterator it = existing_points.find(int_pt); + unsigned int internal_ind; -template -void -DataOutBase::DataOutFilter::write_cell( - unsigned int index, - unsigned int start, - unsigned int d1, - unsigned int d2, - unsigned int d3) -{ - unsigned int base_entry = index * GeometryInfo::vertices_per_cell; - vertices_per_cell = GeometryInfo::vertices_per_cell; - internal_add_cell(base_entry+0, start); - if (dim>=1) - { - internal_add_cell(base_entry+1, start+d1); - if (dim>=2) - { - internal_add_cell(base_entry+2, start+d2+d1); - internal_add_cell(base_entry+3, start+d2); - if (dim>=3) - { - internal_add_cell(base_entry+4, start+d3); - internal_add_cell(base_entry+5, start+d3+d1); - internal_add_cell(base_entry+6, start+d3+d2+d1); - internal_add_cell(base_entry+7, start+d3+d2); - } - } - } -} + // If the point isn't in the set, or we're not filtering duplicate points, add it + if (it == existing_points.end() || !flags.filter_duplicate_vertices) + { + internal_ind = existing_points.size(); + existing_points.insert(std::make_pair(int_pt, internal_ind)); + } + else + { + internal_ind = it->second; + } + // Now add the index to the list of filtered points + filtered_points[index] = internal_ind; + } -void DataOutBase::DataOutFilter::write_data_set(const std::string &name, const unsigned int &dimension, const unsigned int &set_num, const Table<2,double> &data_vectors) -{ - unsigned int num_verts = existing_points.size(); - unsigned int i, r, d, new_dim; - // HDF5/XDMF output only supports 1D or 3D output, so force rearrangement if needed - if (flags.xdmf_hdf5_output && dimension != 1) new_dim = 3; - else new_dim = dimension; - // Record the data set name, dimension, and allocate space for it - data_set_names.push_back(name); - data_set_dims.push_back(new_dim); - data_sets.emplace_back(new_dim*num_verts); + void + DataOutFilter::internal_add_cell(const unsigned int cell_index, + const unsigned int pt_index) + { + filtered_cells[cell_index] = filtered_points[pt_index]; + } - // TODO: averaging, min/max, etc for merged vertices - for (i=0; i &node_data) const + { + node_data.resize(existing_points.size()*node_dim); + + for (Map3DPoint::const_iterator it=existing_points.begin(); it!=existing_points.end(); ++it) + { + for (unsigned int d=0; dsecond+d] = it->first(d); + } + } + + + + void + DataOutFilter::fill_cell_data(const unsigned int local_node_offset, + std::vector &cell_data) const + { + cell_data.resize(filtered_cells.size()); + + for (std::map::const_iterator it=filtered_cells.begin(); it!=filtered_cells.end(); ++it) + { + cell_data[it->first] = it->second+local_node_offset; + } + } + + + + std::string + DataOutFilter::get_data_set_name(const unsigned int set_num) const + { + return data_set_names.at(set_num); + } + + + + unsigned int + DataOutFilter::get_data_set_dim(const unsigned int set_num) const + { + return data_set_dims.at(set_num); + } + + + + const double * + DataOutFilter::get_data_set(const unsigned int set_num) const + { + return &data_sets[set_num][0]; + } + + + + unsigned int + DataOutFilter::n_nodes() const + { + return existing_points.size(); + } + + + + unsigned int + DataOutFilter::n_cells() const + { + return filtered_cells.size()/vertices_per_cell; + } + + + + unsigned int + DataOutFilter::n_data_sets() const + { + return data_set_names.size(); + } + + + + void + DataOutFilter::flush_points () + {} + + + + void + DataOutFilter::flush_cells () + {} + + + + template + void + DataOutFilter::write_cell(const unsigned int index, + const unsigned int start, + const unsigned int d1, + const unsigned int d2, + const unsigned int d3) + { + const unsigned int base_entry = index * GeometryInfo::vertices_per_cell; + vertices_per_cell = GeometryInfo::vertices_per_cell; + internal_add_cell(base_entry+0, start); + if (dim>=1) + { + internal_add_cell(base_entry+1, start+d1); + if (dim>=2) + { + internal_add_cell(base_entry+2, start+d2+d1); + internal_add_cell(base_entry+3, start+d2); + if (dim>=3) + { + internal_add_cell(base_entry+4, start+d3); + internal_add_cell(base_entry+5, start+d3+d1); + internal_add_cell(base_entry+6, start+d3+d2+d1); + internal_add_cell(base_entry+7, start+d3+d2); + } + } + } + } + + + + void DataOutFilter::write_data_set(const std::string &name, + const unsigned int dimension, + const unsigned int set_num, + const Table<2,double> &data_vectors) + { + unsigned int new_dim; + + // HDF5/XDMF output only supports 1D or 3D output, so force rearrangement if needed + if (flags.xdmf_hdf5_output && dimension != 1) + new_dim = 3; + else + new_dim = dimension; + + // Record the data set name, dimension, and allocate space for it + data_set_names.push_back(name); + data_set_dims.push_back(new_dim); + data_sets.emplace_back(new_dim * existing_points.size()); + + // TODO: averaging, min/max, etc for merged vertices + for (unsigned int i=0; i +#include +#include +#include +#include +#include + + + + +template +void +test () +{ + Triangulation tria(Triangulation::MeshSmoothing::none,true); + GridGenerator::hyper_cube(tria, 1., 1. + 1000.0 * std::numeric_limits::epsilon()); + + FE_Q fe1(1); + + DoFHandler dof1(tria); + dof1.distribute_dofs(fe1); + + Vector v1(dof1.n_dofs()); + for (unsigned int i=0; i data_out; + data_out.add_data_vector (dof1, v1, "linear"); + data_out.build_patches (); + + DataOutBase::DataOutFilter data_filter + (DataOutBase::DataOutFilterFlags (true, false)); + + data_out.write_filtered_data (data_filter); + + deallog << "Number of filtered nodes: " << data_filter.n_nodes() << std::endl; + + Triangulation tria2(Triangulation::MeshSmoothing::none,true); + GridGenerator::hyper_cube(tria2, 1., 1 + 10.0 * std::numeric_limits::epsilon()); + + FE_Q fe2(1); + + DoFHandler dof2(tria2); + dof2.distribute_dofs(fe2); + + Vector v2(dof2.n_dofs()); + for (unsigned int i=0; i data_out2; + data_out2.add_data_vector (dof2, v2, "linear"); + data_out2.build_patches (); + + DataOutBase::DataOutFilter data_filter2 + (DataOutBase::DataOutFilterFlags (true, false)); + + data_out2.write_filtered_data (data_filter2); + + deallog << "Number of filtered nodes: " << data_filter2.n_nodes() << std::endl; + + deallog << "ok" << std::endl; +} + + +int main(int argc, char *argv[]) +{ + initlog(); + test<2>(); + test<3>(); + + return 0; +} diff --git a/tests/numerics/data_out_filter_01.with_hdf5=true.output b/tests/numerics/data_out_filter_01.with_hdf5=true.output new file mode 100644 index 0000000000..69858ee394 --- /dev/null +++ b/tests/numerics/data_out_filter_01.with_hdf5=true.output @@ -0,0 +1,7 @@ + +DEAL::Number of filtered nodes: 4 +DEAL::Number of filtered nodes: 1 +DEAL::ok +DEAL::Number of filtered nodes: 8 +DEAL::Number of filtered nodes: 1 +DEAL::ok