From: Wolfgang Bangerth Date: Thu, 31 Aug 2017 21:22:53 +0000 (-0600) Subject: Specialize DataOutBase::Patch for dim==0. X-Git-Tag: v9.0.0-rc1~1133^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=a48d012d7ea4895d05c6be7f0e68c90a95046a50;p=dealii.git Specialize DataOutBase::Patch for dim==0. For dim==1, we output only points, and in this case we can avoid allocating a fair share of the memory that we use in the general case. In reference to #4695. --- diff --git a/include/deal.II/base/data_out_base.h b/include/deal.II/base/data_out_base.h index c165694f15..0aa48c01ba 100644 --- a/include/deal.II/base/data_out_base.h +++ b/include/deal.II/base/data_out_base.h @@ -269,8 +269,8 @@ namespace DataOutBase 1]; /** - * Number of this patch. Since we are not sure patches are handled in the - * same order, always, we better store this. + * Number of this patch. Since we are not sure patches are always + * handled in the same order, we better store this. */ unsigned int patch_index; @@ -347,6 +347,149 @@ namespace DataOutBase * Value to be used if this patch has no neighbor on one side. */ static const unsigned int no_neighbor = numbers::invalid_unsigned_int; + + /** + * @addtogroup Exceptions + * @{ + */ + + /** + * Exception + */ + DeclException2 (ExcInvalidCombinationOfDimensions, + int, int, + << "It is not possible to have a structural dimension of " << arg1 + << " to be larger than the space dimension of the surrounding" + << " space " << arg2); + //@} + }; + + + + /** + * A specialization of the general Patch template that is + * tailored to the case of points, i.e., zero-dimensional objects embedded + * in @p spacedim dimensional space. + * + * The current class is compatible with the general template to allow for + * using the same functions accessing patches of arbitrary dimensionality + * in a generic way. However, it makes some variables that are nonsensical + * for zero-dimensional patches into @p static variables that exist only + * once in the entire program, as opposed to once per patch. Specifically, + * this is the case for the @p neighbors array and the @p n_subdivisions + * member variable that make no sense for zero-dimensional patches because + * points have no natural neighbors across their non-existent faces, nor + * can they reasonably be subdivided. + * + * @author Wolfgang Bangerth, 2017. + */ + template + struct Patch<0,spacedim> + { + /** + * Make the spacedim template parameter available. + */ + static const unsigned int space_dim=spacedim; + + /** + * Corner points of a patch. For the current class of zero-dimensional + * patches, there is of course only a single vertex. + * + * If points_are_available==true, then + * the coordinates of the point at which output is to be generated + * is attached as an additional row to the data table. + */ + Point vertices[1]; + + /** + * An unused, @p static variable that exists only to allow access + * from general code in a generic fashion. + */ + static unsigned int neighbors[1]; + + /** + * Number of this patch. Since we are not sure patches are always + * handled in the same order, we better store this. + */ + unsigned int patch_index; + + /** + * Number of subdivisions with which this patch is to be written. + * 1 means no subdivision, 2 means bisection, 3 + * trisection, etc. + * + * Since subdivision makes no sense for zero-dimensional patches, + * this variable is not used but exists only to allow access + * from general code in a generic fashion. + */ + static unsigned int n_subdivisions; + + /** + * Data vectors. The format is as follows: data(i,.) denotes the + * data belonging to the ith data vector. data.n_cols() + * therefore equals the number of output points; this number is + * of course one for the current class, given that we produce output on + * points. data.n_rows() equals the number of + * data vectors. For the current purpose, a data vector equals one scalar, + * even if multiple scalars may later be interpreted as vectors. + * + * Within each column, data(.,j) are the data values at the + * output point j; for the current class, @p j can only + * be zero. + * + * Since the number of data vectors is usually the same for all patches to + * be printed, data.size() should yield the same value for all + * patches provided. The exception are patches for which + * points_are_available are set, where the actual coordinates of the point + * are appended to the 'data' field, see the documentation of the + * points_are_available flag. + */ + Table<2,float> data; + + /** + * A flag indicating whether the coordinates of the interior patch points + * (assuming that the patch is supposed to be subdivided further) are + * appended to the @p data table (@p true) or not (@p false). The latter + * is the default and in this case the locations of the points interior to + * this patch are computed by (bi-, tri-)linear interpolation from the + * vertices of the patch. + * + * This option exists since patch points may be evaluated using a Mapping + * (rather than by a linear interpolation) and therefore have to be stored + * in the Patch structure. + */ + bool points_are_available; + + /** + * Default constructor. Sets #points_are_available + * to false, and #patch_index to #no_neighbor. + */ + Patch (); + + /** + * Compare the present patch for equality with another one. This is used + * in a few of the automated tests in our testsuite. + */ + bool operator == (const Patch &patch) const; + + /** + * Return an estimate for the memory consumption, in bytes, of this + * object. This is not exact (but will usually be close) because + * calculating the memory usage of trees (e.g., std::map) is + * difficult. + */ + std::size_t memory_consumption () const; + + /** + * Swap the current object's contents with those of the given argument. + */ + void swap (Patch<0,spacedim> &other_patch); + + /** + * Value to be used if this patch has no neighbor on one side. + */ + static const unsigned int no_neighbor = numbers::invalid_unsigned_int; + /** * @addtogroup Exceptions * @{ diff --git a/source/base/data_out_base.cc b/source/base/data_out_base.cc index 3168531d4f..3ec7dbf731 100644 --- a/source/base/data_out_base.cc +++ b/source/base/data_out_base.cc @@ -1533,12 +1533,12 @@ namespace namespace DataOutBase { - template - const unsigned int Patch::space_dim; - const unsigned int Deal_II_IntermediateFlags::format_version = 3; + template + const unsigned int Patch::space_dim; + template const unsigned int Patch::no_neighbor; @@ -1550,10 +1550,8 @@ namespace DataOutBase patch_index(no_neighbor), n_subdivisions (1), points_are_available(false) - // all the other data has a - // constructor of its own, except - // for the "neighbors" field, which - // we set to invalid values. + // all the other data has a constructor of its own, except for the + // "neighbors" field, which we set to invalid values. { for (unsigned int i=0; i::faces_per_cell; ++i) neighbors[i] = no_neighbor; @@ -1638,6 +1636,91 @@ namespace DataOutBase + template + const unsigned int Patch<0,spacedim>::space_dim; + + + template + const unsigned int Patch<0,spacedim>::no_neighbor; + + + template + unsigned int Patch<0,spacedim>::neighbors[1] = { Patch<0,spacedim>::no_neighbor }; + + template + unsigned int Patch<0,spacedim>::n_subdivisions = 1; + + template + Patch<0,spacedim>::Patch () + : + patch_index(no_neighbor), + points_are_available(false) + { + Assert (spacedim<=3, ExcNotImplemented()); + } + + + + template + bool + Patch<0,spacedim>::operator == (const Patch &patch) const + { + const unsigned int dim = 0; + +//TODO: make tolerance relative + const double epsilon=3e-16; + for (unsigned int i=0; i::vertices_per_cell; ++i) + if (vertices[i].distance(patch.vertices[i]) > epsilon) + return false; + + if (patch_index != patch.patch_index) + return false; + + if (points_are_available != patch.points_are_available) + return false; + + if (data.n_rows() != patch.data.n_rows()) + return false; + + if (data.n_cols() != patch.data.n_cols()) + return false; + + for (unsigned int i=0; i + std::size_t + Patch<0,spacedim>::memory_consumption () const + { + return (sizeof(vertices) / sizeof(vertices[0]) * + MemoryConsumption::memory_consumption(vertices[0]) + + + MemoryConsumption::memory_consumption(data) + + + MemoryConsumption::memory_consumption(points_are_available)); + } + + + + template + void + Patch<0,spacedim>::swap (Patch<0,spacedim> &other_patch) + { + std::swap (vertices, other_patch.vertices); + std::swap (patch_index, other_patch.patch_index); + data.swap (other_patch.data); + std::swap (points_are_available, other_patch.points_are_available); + } + + + UcdFlags::UcdFlags (const bool write_preamble) : write_preamble (write_preamble) @@ -7769,8 +7852,7 @@ namespace DataOutBase } - // then read all the data that is - // in this patch + // then read all the data that is in this patch for (unsigned int i=0; i::vertices_per_cell; ++i) in >> patch.vertices[GeometryInfo::ucd_to_deal[i]];