]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Get rid of TriaObject.
authorWolfgang Bangerth <bangerth@colostate.edu>
Tue, 2 Jun 2020 20:35:42 +0000 (14:35 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Tue, 2 Jun 2020 20:35:42 +0000 (14:35 -0600)
TriaObject was a class that was a barely disguised ArrayView, but the class
name wouldn't have told you that. It's not actually very difficult to replace
it by ArrayView, which is what this patch does.

In the process, rename the TriaObjects::get_object() function to
TriaObjects::get_bounding_object_indices() which is both a more
descriptive name of what it returns, and matches the recently
introduced set_bounding_object_indices() function.

include/deal.II/grid/tria_accessor.h
include/deal.II/grid/tria_accessor.templates.h
include/deal.II/grid/tria_faces.h
include/deal.II/grid/tria_levels.h
include/deal.II/grid/tria_object.h [deleted file]
include/deal.II/grid/tria_objects.h
source/grid/tria.cc
source/grid/tria_accessor.cc

index 82998f111bd3bf1ad3f6a2e3035003119c633bc3..862002957e4e7387b65b038a0481841b168eabf2 100644 (file)
@@ -1629,14 +1629,15 @@ private:
    * denotes the indices of the vertices that bound it. And so on.
    */
   void
-  set_bounding_object_indices(const std::initializer_list<int> &o) const;
+  set_bounding_object_indices(
+    const std::initializer_list<int> &new_indices) const;
 
   /**
    * The same as above but for `unsigned int`.
    */
   void
   set_bounding_object_indices(
-    const std::initializer_list<unsigned int> &o) const;
+    const std::initializer_list<unsigned int> &new_indices) const;
 
   /**
    * Set the flag indicating, what <code>line_orientation()</code> will
index bbd52503eef40d3f11958496fc97fcd5ba1a5cf4..2df91d87a17ef33d4a33ca7dbb08ef46718461da 100644 (file)
@@ -525,7 +525,8 @@ namespace internal
       line_index(const TriaAccessor<2, dim, spacedim> &accessor,
                  const unsigned int                    i)
       {
-        return accessor.objects().get_object(accessor.present_index).face(i);
+        return accessor.objects().get_bounding_object_indices(
+          accessor.present_index)[i];
       }
 
 
@@ -570,8 +571,7 @@ namespace internal
                  const unsigned int                    i)
       {
         return accessor.tria->levels[accessor.present_level]
-          ->cells.get_object(accessor.present_index)
-          .face(i);
+          ->cells.get_bounding_object_indices(accessor.present_index)[i];
       }
 
 
@@ -1006,9 +1006,8 @@ namespace internal
       vertex_index(const TriaAccessor<1, dim, spacedim> &accessor,
                    const unsigned int                    corner)
       {
-        return accessor.objects()
-          .get_object(accessor.present_index)
-          .face(corner);
+        return accessor.objects().get_bounding_object_indices(
+          accessor.present_index)[corner];
       }
 
 
index b89e0d54a830d4fc0059f5ea990018cf3494a8ee..7b5445b4e6677c7a7dd9093b86648dfd8d88dbcb 100644 (file)
@@ -18,7 +18,6 @@
 
 #include <deal.II/base/config.h>
 
-#include <deal.II/grid/tria_object.h>
 #include <deal.II/grid/tria_objects.h>
 
 
index 251d3cfbaa64b957e319b0b8e1c2fbc696f999e5..b5ea3b01f3e9ff053cf94755c1f447369e6a85bd 100644 (file)
@@ -21,7 +21,6 @@
 
 #include <deal.II/base/point.h>
 
-#include <deal.II/grid/tria_object.h>
 #include <deal.II/grid/tria_objects.h>
 
 #include <boost/serialization/utility.hpp>
diff --git a/include/deal.II/grid/tria_object.h b/include/deal.II/grid/tria_object.h
deleted file mode 100644 (file)
index c7d6ed5..0000000
+++ /dev/null
@@ -1,112 +0,0 @@
-// ---------------------------------------------------------------------
-//
-// Copyright (C) 1998 - 2020 by the deal.II authors
-//
-// This file is part of the deal.II library.
-//
-// The deal.II library is free software; you can use it, redistribute
-// it, and/or modify it under the terms of the GNU Lesser General
-// Public License as published by the Free Software Foundation; either
-// version 2.1 of the License, or (at your option) any later version.
-// The full text of the license can be found in the file LICENSE.md at
-// the top level directory of deal.II.
-//
-// ---------------------------------------------------------------------
-
-#ifndef dealii_tria_object_h
-#define dealii_tria_object_h
-
-
-#include <deal.II/base/config.h>
-
-#include <deal.II/base/array_view.h>
-#include <deal.II/base/exceptions.h>
-#include <deal.II/base/geometry_info.h>
-
-DEAL_II_NAMESPACE_OPEN
-
-namespace internal
-{
-  namespace TriangulationImplementation
-  {
-    /**
-     * View onto the current geometric object and its faces.
-     *
-     * @note The geometric objects are not saved as separate instances of
-     * TriaObject but in a single vector in TriaObjects.
-     */
-    class TriaObject
-    {
-    public:
-      /**
-       * Constructor.
-       *
-       * @param faces ArrayView inside of the vector in TriaObjects.
-       */
-      TriaObject(const ArrayView<int> &faces)
-        : faces(faces)
-      {}
-
-      /**
-       * Store the content of @p other in the vector of TriaObjects.
-       */
-      TriaObject &
-      operator=(const std::initializer_list<int> &other)
-      {
-        AssertDimension(faces.size(), other.size());
-
-        const std::vector<int> other_v = other;
-
-        for (unsigned int i = 0; i < faces.size(); ++i)
-          faces[i] = other_v[i];
-
-        return *this;
-      }
-
-      /**
-       * The same as above but for `unsigned int`.
-       */
-      TriaObject &
-      operator=(const std::initializer_list<unsigned int> &other)
-      {
-        AssertDimension(faces.size(), other.size());
-
-        const std::vector<unsigned int> other_v = other;
-
-        for (unsigned int i = 0; i < faces.size(); ++i)
-          faces[i] = other_v[i];
-
-        return *this;
-      }
-
-      /**
-       * Return index of the @p i-th face.
-       */
-      int
-      face(const unsigned int i) const
-      {
-        return faces[i];
-      }
-
-      /**
-       * Set index of the @p i-th face.
-       */
-      void
-      set_face(const unsigned int i, const int index)
-      {
-        faces[i] = index;
-      }
-
-    private:
-      /**
-       * List of face this object is made up.
-       */
-      const ArrayView<int> faces;
-    };
-  } // namespace TriangulationImplementation
-} // namespace internal
-
-
-DEAL_II_NAMESPACE_CLOSE
-
-#endif
index 592224a442dc93bbc4cde210be7f05f1580fa7ab..982bd75a862e9491287544e767ef43e8f4fdd3b6 100644 (file)
@@ -22,8 +22,6 @@
 #include <deal.II/base/exceptions.h>
 #include <deal.II/base/geometry_info.h>
 
-#include <deal.II/grid/tria_object.h>
-
 #include <vector>
 
 DEAL_II_NAMESPACE_OPEN
@@ -82,10 +80,14 @@ namespace internal
       n_objects() const;
 
       /**
-       * Return a view on the @p index-th geometric object.
+       * Return a view on the indices of the objects that bound the @p
+       * index-th object stored by the current object. For example, if
+       * the current object stores cells, then this function returns
+       * the equivalent of an array containing the indices of the
+       * faces that bound the @p index-th cell.
        */
-      TriaObject
-      get_object(const unsigned int index);
+      ArrayView<int>
+      get_bounding_object_indices(const unsigned int index);
 
       /**
        * Index of the even children of an object. Since when objects are
@@ -389,8 +391,8 @@ namespace internal
 
 
 
-    inline TriaObject
-    TriaObjects::get_object(const unsigned int index)
+    inline ArrayView<int>
+    TriaObjects::get_bounding_object_indices(const unsigned int index)
     {
       // assume that each cell has the same number of faces
 
@@ -405,8 +407,8 @@ namespace internal
       else
         AssertThrow(false, ExcNotImplemented());
 
-      return {
-        ArrayView<int>(cells.data() + index * faces_per_cell, faces_per_cell)};
+      return ArrayView<int>(cells.data() + index * faces_per_cell,
+                            faces_per_cell);
     }
 
 
index 7ee5647dba93b7d80ed9f483037c5200cc909681..4540facbfd821109c7cb5e73ae2853f686fa642c 100644 (file)
@@ -4365,17 +4365,18 @@ namespace internal
                                   for (const unsigned int q :
                                        GeometryInfo<dim>::face_indices())
                                     {
-                                      const int index = triangulation.levels[l]
-                                                          ->cells.get_object(h)
-                                                          .face(q);
+                                      const int index =
+                                        triangulation.levels[l]
+                                          ->cells.get_bounding_object_indices(
+                                            h)[q];
                                       if (index == switch_1_index)
                                         triangulation.levels[l]
-                                          ->cells.get_object(h)
-                                          .set_face(q, switch_2_index);
+                                          ->cells.get_bounding_object_indices(
+                                            h)[q] = switch_2_index;
                                       else if (index == switch_2_index)
                                         triangulation.levels[l]
-                                          ->cells.get_object(h)
-                                          .set_face(q, switch_1_index);
+                                          ->cells.get_bounding_object_indices(
+                                            h)[q] = switch_1_index;
                                     }
                               // now we have to copy
                               // all information of the
@@ -6338,14 +6339,16 @@ namespace internal
                                    ++l)
                                 {
                                   const int this_index =
-                                    triangulation.faces->quads.get_object(q)
-                                      .face(l);
+                                    triangulation.faces->quads
+                                      .get_bounding_object_indices(q)[l];
                                   if (this_index == old_index_0)
-                                    triangulation.faces->quads.get_object(q)
-                                      .set_face(l, new_index_0);
+                                    triangulation.faces->quads
+                                      .get_bounding_object_indices(q)[l] =
+                                      new_index_0;
                                   else if (this_index == old_index_1)
-                                    triangulation.faces->quads.get_object(q)
-                                      .set_face(l, new_index_1);
+                                    triangulation.faces->quads
+                                      .get_bounding_object_indices(q)[l] =
+                                      new_index_1;
                                 }
                             // now we have to copy all information of
                             // the two lines
@@ -6437,16 +6440,16 @@ namespace internal
                                   {
                                     const int face_index =
                                       triangulation.levels[l]
-                                        ->cells.get_object(h)
-                                        .face(q);
+                                        ->cells.get_bounding_object_indices(
+                                          h)[q];
                                     if (face_index == switch_1_index)
                                       triangulation.levels[l]
-                                        ->cells.get_object(h)
-                                        .set_face(q, switch_2_index);
+                                        ->cells.get_bounding_object_indices(
+                                          h)[q] = switch_2_index;
                                     else if (face_index == switch_2_index)
                                       triangulation.levels[l]
-                                        ->cells.get_object(h)
-                                        .set_face(q, switch_1_index);
+                                        ->cells.get_bounding_object_indices(
+                                          h)[q] = switch_1_index;
                                   }
                             // now we have to copy all information of
                             // the two quads
index 8fadc244800d3e158513b54248cbb876ec8fd7fd..51a0971543368685c6536be080fad4f1b5877229 100644 (file)
@@ -1489,9 +1489,19 @@ const unsigned int
 template <int structdim, int dim, int spacedim>
 void
 TriaAccessor<structdim, dim, spacedim>::set_bounding_object_indices(
-  const std::initializer_list<int> &object) const
+  const std::initializer_list<int> &new_indices) const
 {
-  this->objects().get_object(this->present_index) = object;
+  const ArrayView<int> bounding_object_index_ref =
+    this->objects().get_bounding_object_indices(this->present_index);
+
+  AssertDimension(bounding_object_index_ref.size(), new_indices.size());
+
+  unsigned int i = 0;
+  for (const auto &new_index : new_indices)
+    {
+      bounding_object_index_ref[i] = new_index;
+      ++i;
+    }
 }
 
 
@@ -1499,9 +1509,19 @@ TriaAccessor<structdim, dim, spacedim>::set_bounding_object_indices(
 template <int structdim, int dim, int spacedim>
 void
 TriaAccessor<structdim, dim, spacedim>::set_bounding_object_indices(
-  const std::initializer_list<unsigned int> &object) const
+  const std::initializer_list<unsigned int> &new_indices) const
 {
-  this->objects().get_object(this->present_index) = object;
+  const ArrayView<int> bounding_object_index_ref =
+    this->objects().get_bounding_object_indices(this->present_index);
+
+  AssertDimension(bounding_object_index_ref.size(), new_indices.size());
+
+  unsigned int i = 0;
+  for (const auto &new_index : new_indices)
+    {
+      bounding_object_index_ref[i] = new_index;
+      ++i;
+    }
 }
 
 

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.