]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add ReferenceCell::max_n_vertices<dim>(). 17974/head
authorDavid Wells <drwells@email.unc.edu>
Sun, 5 Jan 2025 18:10:20 +0000 (13:10 -0500)
committerDavid Wells <drwells@email.unc.edu>
Sun, 5 Jan 2025 20:29:23 +0000 (15:29 -0500)
Like the other max functions, this helps us semantically disambiguate
between a maximum value and the hypercube value.

I had to add more preprocessor checks to work around Windows Server
2019's MSVC's problems with boost. I'm not sure when GitHub will start
removing support for that OS, but it may be soon: it entered extended
support a year ago [1]. I think we should just leave these preprocessor
defines until we retire support for that OS, which may be after the next
release.

[1] https://learn.microsoft.com/en-us/lifecycle/products/windows-server-2019

14 files changed:
include/deal.II/fe/mapping.h
include/deal.II/fe/mapping_fe_field.h
include/deal.II/fe/mapping_q1_eulerian.h
include/deal.II/fe/mapping_q_cache.h
include/deal.II/fe/mapping_q_eulerian.h
include/deal.II/grid/reference_cell.h
include/deal.II/grid/tria_accessor.h
source/fe/mapping.cc
source/fe/mapping_fe_field.cc
source/fe/mapping_q.cc
source/fe/mapping_q1_eulerian.cc
source/fe/mapping_q_cache.cc
source/fe/mapping_q_eulerian.cc
source/grid/tria.cc

index e578c853547e264f1c2c95e7015c51e48fc29466..7681b7765d55feb5c7718263d18b5522d4e693b4 100644 (file)
@@ -349,7 +349,12 @@ public:
    * <code>cell-@>vertex(v)</code>.
    */
   virtual boost::container::small_vector<Point<spacedim>,
-                                         GeometryInfo<dim>::vertices_per_cell>
+#ifndef _MSC_VER
+                                         ReferenceCell::max_n_vertices<dim>()
+#else
+                                         GeometryInfo<dim>::vertices_per_cell
+#endif
+                                         >
   get_vertices(
     const typename Triangulation<dim, spacedim>::cell_iterator &cell) const;
 
@@ -362,7 +367,12 @@ public:
    * @param[in] face_no The number of the face within the cell.
    */
   boost::container::small_vector<Point<spacedim>,
-                                 GeometryInfo<dim>::vertices_per_face>
+#ifndef _MSC_VER
+                                 ReferenceCell::max_n_vertices<dim - 1>()
+#else
+                                 GeometryInfo<dim - 1>::vertices_per_cell
+#endif
+                                 >
   get_vertices(const typename Triangulation<dim, spacedim>::cell_iterator &cell,
                const unsigned int face_no) const;
 
index 574999ee659c305085c965091ec2e9425bc929a4..1926de6bb7786961c223eeeddb067f7f90355fda 100644 (file)
@@ -170,7 +170,12 @@ public:
    * that was passed at construction time.
    */
   virtual boost::container::small_vector<Point<spacedim>,
-                                         GeometryInfo<dim>::vertices_per_cell>
+#ifndef _MSC_VER
+                                         ReferenceCell::max_n_vertices<dim>()
+#else
+                                         GeometryInfo<dim>::vertices_per_cell
+#endif
+                                         >
   get_vertices(const typename Triangulation<dim, spacedim>::cell_iterator &cell)
     const override;
 
index 0bf4e6a9f956156a4b7d75188687ad0f7f18e7b6..7810f86ee98233584576501f8214c58ab6946a4f 100644 (file)
@@ -118,7 +118,12 @@ public:
    * addition to the geometry of the cell.
    */
   virtual boost::container::small_vector<Point<spacedim>,
-                                         GeometryInfo<dim>::vertices_per_cell>
+#ifndef _MSC_VER
+                                         ReferenceCell::max_n_vertices<dim>()
+#else
+                                         GeometryInfo<dim>::vertices_per_cell
+#endif
+                                         >
   get_vertices(const typename Triangulation<dim, spacedim>::cell_iterator &cell)
     const override;
 
index 99dfaf3020dd31dd6a2794fb92990735f3b564f5..a7a6398af79a8c39384315f5c591d01f806df62f 100644 (file)
@@ -195,7 +195,12 @@ public:
    * @copydoc Mapping::get_vertices()
    */
   virtual boost::container::small_vector<Point<spacedim>,
-                                         GeometryInfo<dim>::vertices_per_cell>
+#ifndef _MSC_VER
+                                         ReferenceCell::max_n_vertices<dim>()
+#else
+                                         GeometryInfo<dim>::vertices_per_cell
+#endif
+                                         >
   get_vertices(const typename Triangulation<dim, spacedim>::cell_iterator &cell)
     const override;
 
index 4ca9953b03e0aee93d746c5f0016452aa7549fe2..f6355b86736f22890307290a6ed7b368f785faea 100644 (file)
@@ -122,7 +122,12 @@ public:
    * addition to the geometry of the cell.
    */
   virtual boost::container::small_vector<Point<spacedim>,
-                                         GeometryInfo<dim>::vertices_per_cell>
+#ifndef _MSC_VER
+                                         ReferenceCell::max_n_vertices<dim>()
+#else
+                                         GeometryInfo<dim>::vertices_per_cell
+#endif
+                                         >
   get_vertices(const typename Triangulation<dim, spacedim>::cell_iterator &cell)
     const override;
 
index 64268f40a6572880caa83d92f51e1af048006044..98b363a0cf67b99c719ebbe0aa6ab5fce5b42c7b 100644 (file)
@@ -285,6 +285,17 @@ public:
   unsigned int
   n_vertices() const;
 
+  /**
+   * Return the maximum number of vertices an object of dimension `structdim`
+   * can have. This is always the number of vertices of a
+   * `structdim`-dimensional hypercube.
+   *
+   * @see ReferenceCell::max_n_faces()
+   */
+  template <int structdim>
+  static constexpr unsigned int
+  max_n_vertices();
+
   /**
    * Return an object that can be thought of as an array containing all
    * indices from zero to n_vertices().
@@ -1591,6 +1602,15 @@ ReferenceCell::n_vertices() const
 
 
 
+template <int structdim>
+inline constexpr unsigned int
+ReferenceCell::max_n_vertices()
+{
+  return GeometryInfo<structdim>::vertices_per_cell;
+}
+
+
+
 inline unsigned int
 ReferenceCell::n_lines() const
 {
index b3ea594ba41fcb6120da45e8bfa96909c85f576e..6807698128a3d638d74307cecb003989b196ab23 100644 (file)
@@ -5311,9 +5311,8 @@ TriaAccessor<structdim, dim, spacedim>::vertex_index(
     {
       // This branch should only be used after the cell vertex index cache is
       // set up
-      constexpr unsigned int max_vertices_per_cell = 1 << dim;
-      const std::size_t      my_index =
-        static_cast<std::size_t>(this->present_index) * max_vertices_per_cell;
+      const auto my_index = static_cast<std::size_t>(this->present_index) *
+                            ReferenceCell::max_n_vertices<dim>();
       AssertIndexRange(my_index + corner,
                        this->tria->levels[this->present_level]
                          ->cell_vertex_indices_cache.size());
@@ -6347,7 +6346,12 @@ double
 TriaAccessor<structdim, dim, spacedim>::diameter() const
 {
   boost::container::small_vector<Point<spacedim>,
-                                 GeometryInfo<structdim>::vertices_per_cell>
+#  ifndef _MSC_VER
+                                 ReferenceCell::max_n_vertices<structdim>()
+#  else
+                                 GeometryInfo<structdim>::vertices_per_cell
+#  endif
+                                 >
     vertices(this->n_vertices());
 
   for (unsigned int v = 0; v < vertices.size(); ++v)
index c7fb0891e249cc456e42baab7cda7dd23491098f..a6386f89d456fb192c862c1dd4fb7e557b42d139 100644 (file)
@@ -34,12 +34,22 @@ DEAL_II_NAMESPACE_OPEN
 
 template <int dim, int spacedim>
 boost::container::small_vector<Point<spacedim>,
-                               GeometryInfo<dim>::vertices_per_cell>
+#  ifndef _MSC_VER
+                               ReferenceCell::max_n_vertices<dim>()
+#  else
+                               GeometryInfo<dim>::vertices_per_cell
+#  endif
+                               >
 Mapping<dim, spacedim>::get_vertices(
   const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
 {
   boost::container::small_vector<Point<spacedim>,
-                                 GeometryInfo<dim>::vertices_per_cell>
+#  ifndef _MSC_VER
+                                 ReferenceCell::max_n_vertices<dim>()
+#  else
+                                 GeometryInfo<dim>::vertices_per_cell
+#  endif
+                                 >
     vertices;
   for (const unsigned int i : cell->vertex_indices())
     vertices.push_back(cell->vertex(i));
@@ -51,13 +61,23 @@ Mapping<dim, spacedim>::get_vertices(
 
 template <int dim, int spacedim>
 boost::container::small_vector<Point<spacedim>,
-                               GeometryInfo<dim>::vertices_per_face>
+#  ifndef _MSC_VER
+                               ReferenceCell::max_n_vertices<dim - 1>()
+#  else
+                               GeometryInfo<dim - 1>::vertices_per_cell
+#  endif
+                               >
 Mapping<dim, spacedim>::get_vertices(
   const typename Triangulation<dim, spacedim>::cell_iterator &cell,
   const unsigned int                                          face_no) const
 {
   boost::container::small_vector<Point<spacedim>,
-                                 GeometryInfo<dim>::vertices_per_face>
+#  ifndef _MSC_VER
+                                 ReferenceCell::max_n_vertices<dim - 1>()
+#  else
+                                 GeometryInfo<dim - 1>::vertices_per_cell
+#  endif
+                                 >
     face_vertices;
 
   const auto &cell_vertices    = get_vertices(cell);
index 4711207f28cf904ce85e16493bced58cb0c6bc55..58d9f8f0cbac6af8bfa1e1bdd8fe7385e5e7e529 100644 (file)
@@ -447,7 +447,12 @@ MappingFEField<dim, spacedim, VectorType>::is_compatible_with(
 
 template <int dim, int spacedim, typename VectorType>
 boost::container::small_vector<Point<spacedim>,
-                               GeometryInfo<dim>::vertices_per_cell>
+#ifndef _MSC_VER
+                               ReferenceCell::max_n_vertices<dim>()
+#else
+                               GeometryInfo<dim>::vertices_per_cell
+#endif
+                               >
 MappingFEField<dim, spacedim, VectorType>::get_vertices(
   const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
 {
@@ -485,7 +490,12 @@ MappingFEField<dim, spacedim, VectorType>::get_vertices(
     uses_level_dofs ? *euler_vector[cell->level()] : *euler_vector[0];
 
   boost::container::small_vector<Point<spacedim>,
-                                 GeometryInfo<dim>::vertices_per_cell>
+#ifndef _MSC_VER
+                                 ReferenceCell::max_n_vertices<dim>()
+#else
+                                 GeometryInfo<dim>::vertices_per_cell
+#endif
+                                 >
     vertices(cell->n_vertices());
   for (unsigned int i = 0; i < dofs_per_cell; ++i)
     {
index 734fb98646513b288fccabc479efb0c5e4f54d6d..55c2fe6aee9df29bf9041740f95a5a89ca6896a1 100644 (file)
@@ -642,7 +642,12 @@ MappingQ<dim, spacedim>::transform_points_real_to_unit_cell(
   AssertDimension(real_points.size(), unit_points.size());
   std::vector<Point<spacedim>> support_points_higher_order;
   boost::container::small_vector<Point<spacedim>,
-                                 GeometryInfo<dim>::vertices_per_cell>
+#ifndef _MSC_VER
+                                 ReferenceCell::max_n_vertices<dim>()
+#else
+                                 GeometryInfo<dim>::vertices_per_cell
+#endif
+                                 >
     vertices;
   if (polynomial_degree == 1)
     vertices = this->get_vertices(cell);
index 0bd74b6c7d32bcc721957d73cef19bfed1cc4f91..79358a1ccb9479efa1dc9ac8da62687e5cd9f205 100644 (file)
@@ -49,13 +49,23 @@ MappingQ1Eulerian<dim, VectorType, spacedim>::MappingQ1Eulerian(
 
 template <int dim, typename VectorType, int spacedim>
 boost::container::small_vector<Point<spacedim>,
-                               GeometryInfo<dim>::vertices_per_cell>
+#ifndef _MSC_VER
+                               ReferenceCell::max_n_vertices<dim>()
+#else
+                               GeometryInfo<dim>::vertices_per_cell
+#endif
+                               >
 MappingQ1Eulerian<dim, VectorType, spacedim>::get_vertices(
   const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
 {
   boost::container::small_vector<Point<spacedim>,
-                                 GeometryInfo<dim>::vertices_per_cell>
-    vertices(GeometryInfo<dim>::vertices_per_cell);
+#ifndef _MSC_VER
+                                 ReferenceCell::max_n_vertices<dim>()
+#else
+                                 GeometryInfo<dim>::vertices_per_cell
+#endif
+                                 >
+    vertices(cell->n_vertices());
   // The assertions can not be in the constructor, since this would
   // require to call dof_handler.distribute_dofs(fe) *before* the mapping
   // object is constructed, which is not necessarily what we want.
index 750a9602d12c2cfff61fa912ac8dd2895668274b..52f0bca873d8acd485dd16d314e9306004404e18 100644 (file)
@@ -728,7 +728,12 @@ MappingQCache<dim, spacedim>::compute_mapping_support_points(
 
 template <int dim, int spacedim>
 boost::container::small_vector<Point<spacedim>,
-                               GeometryInfo<dim>::vertices_per_cell>
+#ifndef _MSC_VER
+                               ReferenceCell::max_n_vertices<dim>()
+#else
+                               GeometryInfo<dim>::vertices_per_cell
+#endif
+                               >
 MappingQCache<dim, spacedim>::get_vertices(
   const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
 {
@@ -742,8 +747,12 @@ MappingQCache<dim, spacedim>::get_vertices(
   AssertIndexRange(cell->index(), (*support_point_cache)[cell->level()].size());
   const auto ptr = (*support_point_cache)[cell->level()][cell->index()].begin();
   return boost::container::small_vector<Point<spacedim>,
-                                        GeometryInfo<dim>::vertices_per_cell>(
-    ptr, ptr + cell->n_vertices());
+#ifndef _MSC_VER
+                                        ReferenceCell::max_n_vertices<dim>()
+#else
+                                        GeometryInfo<dim>::vertices_per_cell
+#endif
+                                        >(ptr, ptr + cell->n_vertices());
 }
 
 
index 63267c1f26caacd6198bb9ce14411a6a1155e5ac..3569730fd151febd0081a3b2930b6068a58be8f1 100644 (file)
@@ -99,7 +99,12 @@ MappingQEulerian<dim, VectorType, spacedim>::SupportQuadrature::
 
 template <int dim, typename VectorType, int spacedim>
 boost::container::small_vector<Point<spacedim>,
-                               GeometryInfo<dim>::vertices_per_cell>
+#ifndef _MSC_VER
+                               ReferenceCell::max_n_vertices<dim>()
+#else
+                               GeometryInfo<dim>::vertices_per_cell
+#endif
+                               >
 MappingQEulerian<dim, VectorType, spacedim>::get_vertices(
   const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
 {
@@ -107,9 +112,13 @@ MappingQEulerian<dim, VectorType, spacedim>::get_vertices(
   const std::vector<Point<spacedim>> a = compute_mapping_support_points(cell);
 
   boost::container::small_vector<Point<spacedim>,
-                                 GeometryInfo<dim>::vertices_per_cell>
-    vertex_locations(a.begin(),
-                     a.begin() + GeometryInfo<dim>::vertices_per_cell);
+#ifndef _MSC_VER
+                                 ReferenceCell::max_n_vertices<dim>()
+#else
+                                 GeometryInfo<dim>::vertices_per_cell
+#endif
+                                 >
+    vertex_locations(a.begin(), a.begin() + cell->n_vertices());
 
   return vertex_locations;
 }
index 90f9b684e3ac0975ed1271941392b740019355b7..109da0ad42032bc7573a430ef1d7f5713a5a4995 100644 (file)
@@ -3498,7 +3498,7 @@ namespace internal
         [[maybe_unused]] unsigned int counter = 0;
 
         std::vector<unsigned int> key;
-        key.reserve(GeometryInfo<structdim>::vertices_per_cell);
+        key.reserve(ReferenceCell::max_n_vertices<structdim>());
 
         for (unsigned int o = 0; o < obj.n_objects(); ++o)
           {
@@ -15953,14 +15953,15 @@ void Triangulation<dim, spacedim>::reset_cell_vertex_indices_cache()
 {
   for (unsigned int l = 0; l < levels.size(); ++l)
     {
-      constexpr unsigned int     max_vertices_per_cell = 1 << dim;
       std::vector<unsigned int> &cache = levels[l]->cell_vertex_indices_cache;
       cache.clear();
-      cache.resize(levels[l]->refine_flags.size() * max_vertices_per_cell,
+      cache.resize(levels[l]->refine_flags.size() *
+                     ReferenceCell::max_n_vertices<dim>(),
                    numbers::invalid_unsigned_int);
       for (const auto &cell : cell_iterators_on_level(l))
         {
-          const unsigned int my_index = cell->index() * max_vertices_per_cell;
+          const unsigned int my_index =
+            cell->index() * ReferenceCell::max_n_vertices<dim>();
 
           // to reduce the cost of this function when passing down into quads,
           // then lines, then vertices, we use a more low-level access method

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.