tria_objects.h

Go to the documentation of this file.

// -----------------------------------------------------------------------------
//
// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception OR LGPL-2.1-or-later
// Copyright (C) 2006 - 2026 by the deal.II authors
//
// This file is part of the deal.II library.
//
// Detailed license information governing the source code and contributions
// can be found in LICENSE.md and CONTRIBUTING.md at the top level directory.
//
// -----------------------------------------------------------------------------
 
#ifndef dealii_tria_objects_h
#define dealii_tria_objects_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>
 
#include <vector>
 
DEAL_II_NAMESPACE_OPEN
 
// Forward declarations
#ifndef DOXYGEN
template <int dim, int spacedim>
DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim<dim, spacedim>))
class Triangulation;
template <class Accessor>
class TriaRawIterator;
template <int, int, int>
class TriaAccessor;
#endif
 
namespace internal
{
  namespace TriangulationImplementation
  {
    class TriaObjects
    {
    public:
      TriaObjects();
 
      TriaObjects(const unsigned int structdim,
                  const unsigned int children_per_object,
                  const unsigned int faces_per_object);
 
      void
      allocate(const std::size_t n_objects);
 
      void
      allocate_end(const std::size_t new_objects_in_pairs,
                   const std::size_t new_objects_single);
 
      unsigned int structdim;
 
      unsigned int children_per_object;
 
      unsigned int faces_per_object;
 
      std::size_t
      n_objects() const;
 
      int
      get_bounding_object(const int index, const unsigned int face_no) const;
 
      void
      set_bounding_object(const int          index,
                          const unsigned int face_no,
                          const int          neighbor_index);
 
      ArrayView<int>
      get_bounding_object_indices(const unsigned int index);
 
      std::vector<int> children;
 
      std::vector<std::uint8_t> refinement_cases;
 
      std::vector<bool> used;
 
      std::vector<bool> user_flags;
 
 
      struct BoundaryOrMaterialId
      {
        union
        {
          types::boundary_id boundary_id;
          types::material_id material_id;
        };
 
 
        BoundaryOrMaterialId();
 
        static std::size_t
        memory_consumption();
 
        template <class Archive>
        void
        serialize(Archive &ar, const unsigned int version);
      };
 
      std::vector<BoundaryOrMaterialId> boundary_or_material_id;
 
      std::vector<types::manifold_id> manifold_id;
 
      template <int structdim, int dim, int spacedim>
      ::TriaRawIterator<::TriaAccessor<structdim, dim, spacedim>>
      next_free_single_object(const Triangulation<dim, spacedim> &tria);
 
      template <int structdim, int dim, int spacedim>
      ::TriaRawIterator<::TriaAccessor<structdim, dim, spacedim>>
      next_free_pair_object(const Triangulation<dim, spacedim> &tria);
 
      template <int dim, int spacedim>
      typename Triangulation<dim, spacedim>::raw_hex_iterator
      next_free_hex(const Triangulation<dim, spacedim> &tria,
                    const unsigned int                  level);
 
      void *&
      user_pointer(const unsigned int i);
 
      const void *
      user_pointer(const unsigned int i) const;
 
      unsigned int &
      user_index(const unsigned int i);
 
      unsigned int
      user_index(const unsigned int i) const;
 
      void
      clear_user_data(const unsigned int i);
 
      void
      clear_user_data();
 
      void
      clear_user_flags();
 
      std::size_t
      memory_consumption() const;
 
      template <class Archive>
      void
      serialize(Archive &ar, const unsigned int version);
 
      DeclException0(ExcPointerIndexClash);
 
      unsigned int next_free_single;
 
      unsigned int next_free_pair;
 
      bool reverse_order_next_free_single;
 
      struct UserData
      {
        union
        {
          void *p;
          unsigned int i;
        };
 
        UserData()
        {
          p = nullptr;
        }
 
        template <class Archive>
        void
        serialize(Archive &ar, const unsigned int version);
      };
 
      enum UserDataType
      {
        data_unknown,
        data_pointer,
        data_index
      };
 
 
      std::vector<UserData> user_data;
 
      mutable UserDataType user_data_type;
 
    private:
      std::vector<int> bounding_objects;
    };
 
 
    //----------------------------------------------------------------------//
 
    inline std::size_t
    TriaObjects::n_objects() const
    {
      // ensure that sizes are consistent, and then return one that
      // corresponds to the number of objects
      AssertDimension(bounding_objects.size(),
                      manifold_id.size() * faces_per_object);
      return manifold_id.size();
    }
 
 
 
    inline int
    TriaObjects::get_bounding_object(const int          index,
                                     const unsigned int face_no) const
    {
      // Accessors with an index of -1 shouldn't get here
      Assert(index != -1, ExcInternalError());
      AssertIndexRange(face_no, faces_per_object);
      const auto i = index * faces_per_object + face_no;
      AssertIndexRange(i, bounding_objects.size());
      return bounding_objects[i];
    }
 
 
 
    inline void
    TriaObjects::set_bounding_object(const int          index,
                                     const unsigned int face_no,
                                     const int          neighbor_index)
    {
      // Accessors with an index of -1 shouldn't get here
      Assert(index != -1, ExcInternalError());
      Assert(neighbor_index >= -1, ExcInternalError());
      AssertIndexRange(face_no, faces_per_object);
      const auto i = index * faces_per_object + face_no;
      AssertIndexRange(i, bounding_objects.size());
      bounding_objects[i] = neighbor_index;
    }
 
 
 
    inline ArrayView<int>
    TriaObjects::get_bounding_object_indices(const unsigned int index)
    {
      // each cell has the same number of subcell objects
      return ArrayView<int>(bounding_objects.data() + index * faces_per_object,
                            faces_per_object);
    }
 
 
 
    inline TriaObjects::BoundaryOrMaterialId::BoundaryOrMaterialId()
    {
      material_id = numbers::invalid_material_id;
    }
 
 
 
    inline std::size_t
    TriaObjects::BoundaryOrMaterialId::memory_consumption()
    {
      return sizeof(BoundaryOrMaterialId);
    }
 
 
 
    template <class Archive>
    void
    TriaObjects::BoundaryOrMaterialId::serialize(Archive &ar,
                                                 const unsigned int /*version*/)
    {
      // serialize this
      // structure by
      // writing and
      // reading the larger
      // of the two values,
      // in order to make
      // sure we get all
      // bits
      if (sizeof(material_id) > sizeof(boundary_id))
        ar &material_id;
      else
        ar &boundary_id;
    }
 
 
    inline void *&
    TriaObjects::user_pointer(const unsigned int i)
    {
      Assert(user_data_type == data_unknown || user_data_type == data_pointer,
             ExcPointerIndexClash());
      user_data_type = data_pointer;
 
      AssertIndexRange(i, user_data.size());
      return user_data[i].p;
    }
 
 
    inline const void *
    TriaObjects::user_pointer(const unsigned int i) const
    {
      Assert(user_data_type == data_unknown || user_data_type == data_pointer,
             ExcPointerIndexClash());
      user_data_type = data_pointer;
 
      AssertIndexRange(i, user_data.size());
      return user_data[i].p;
    }
 
 
    inline unsigned int &
    TriaObjects::user_index(const unsigned int i)
    {
      Assert(user_data_type == data_unknown || user_data_type == data_index,
             ExcPointerIndexClash());
      user_data_type = data_index;
 
      AssertIndexRange(i, user_data.size());
      return user_data[i].i;
    }
 
 
    inline void
    TriaObjects::clear_user_data(const unsigned int i)
    {
      AssertIndexRange(i, user_data.size());
      user_data[i].i = 0;
    }
 
 
    inline TriaObjects::TriaObjects()
      : structdim(numbers::invalid_unsigned_int)
      , children_per_object(numbers::invalid_unsigned_int)
      , faces_per_object(numbers::invalid_unsigned_int)
      , next_free_single(numbers::invalid_unsigned_int)
      , next_free_pair(numbers::invalid_unsigned_int)
      , reverse_order_next_free_single(false)
      , user_data_type(data_unknown)
    {}
 
 
    inline TriaObjects::TriaObjects(const unsigned int structdim,
                                    const unsigned int children_per_object,
                                    const unsigned int faces_per_object)
      : structdim(structdim)
      , children_per_object(children_per_object)
      , faces_per_object(faces_per_object)
      , next_free_single(numbers::invalid_unsigned_int)
      , next_free_pair(numbers::invalid_unsigned_int)
      , reverse_order_next_free_single(false)
      , user_data_type(data_unknown)
    {}
 
 
    inline unsigned int
    TriaObjects::user_index(const unsigned int i) const
    {
      Assert(user_data_type == data_unknown || user_data_type == data_index,
             ExcPointerIndexClash());
      user_data_type = data_index;
 
      AssertIndexRange(i, user_data.size());
      return user_data[i].i;
    }
 
 
    inline void
    TriaObjects::clear_user_data()
    {
      user_data_type = data_unknown;
      for (auto &data : user_data)
        data.p = nullptr;
    }
 
 
    inline void
    TriaObjects::clear_user_flags()
    {
      user_flags.assign(user_flags.size(), false);
    }
 
 
    template <class Archive>
    void
    TriaObjects::UserData::serialize(Archive &ar, const unsigned int)
    {
      // serialize this as an integer
      ar &i;
    }
 
 
 
    template <class Archive>
    void
    TriaObjects::serialize(Archive &ar, const unsigned int)
    {
      ar                                   &structdim;
      ar                                   &children_per_object;
      ar                                   &faces_per_object;
      ar &bounding_objects                 &children;
      ar                                   &refinement_cases;
      ar                                   &used;
      ar                                   &user_flags;
      ar                                   &boundary_or_material_id;
      ar                                   &manifold_id;
      ar &next_free_single &next_free_pair &reverse_order_next_free_single;
      ar &user_data                        &user_data_type;
    }
 
 
    //----------------------------------------------------------------------//
 
    template <int structdim_, int dim, int spacedim>
    ::TriaRawIterator<::TriaAccessor<structdim_, dim, spacedim>>
    TriaObjects::next_free_single_object(
      const Triangulation<dim, spacedim> &tria)
    {
      // TODO: Think of a way to ensure that we are using the correct
      // triangulation, i.e. the one containing *this.
 
      AssertDimension(structdim_, this->structdim);
 
      int pos = next_free_single, last = used.size() - 1;
      if (!reverse_order_next_free_single)
        {
          // first sweep forward, only use really single slots, do not use
          // pair slots
          for (; pos < last; ++pos)
            if (!used[pos])
              if (used[++pos])
                {
                  // this was a single slot
                  pos -= 1;
                  break;
                }
          if (pos >= last)
            {
              reverse_order_next_free_single = true;
              next_free_single               = used.size() - 1;
              pos                            = used.size() - 1;
            }
          else
            next_free_single = pos + 1;
        }
 
      if (reverse_order_next_free_single)
        {
          // second sweep, use all slots, even
          // in pairs
          for (; pos >= 0; --pos)
            if (!used[pos])
              break;
          if (pos > 0)
            next_free_single = pos - 1;
          else
            // no valid single object anymore
            return ::TriaRawIterator<
              ::TriaAccessor<structdim_, dim, spacedim>>(&tria, -1, -1);
        }
 
      return ::TriaRawIterator<
        ::TriaAccessor<structdim_, dim, spacedim>>(&tria, 0, pos);
    }
 
 
 
    template <int structdim_, int dim, int spacedim>
    ::TriaRawIterator<::TriaAccessor<structdim_, dim, spacedim>>
    TriaObjects::next_free_pair_object(const Triangulation<dim, spacedim> &tria)
    {
      // TODO: Think of a way to ensure that we are using the correct
      // triangulation, i.e. the one containing *this.
 
      AssertDimension(structdim_, this->structdim);
 
      int pos = next_free_pair, last = used.size() - 1;
      for (; pos < last; ++pos)
        if (!used[pos])
          if (!used[++pos])
            {
              // this was a pair slot
              pos -= 1;
              break;
            }
      if (pos >= last)
        // no free slot
        return ::TriaRawIterator<
          ::TriaAccessor<structdim_, dim, spacedim>>(&tria, -1, -1);
      else
        next_free_pair = pos + 2;
 
      return ::TriaRawIterator<
        ::TriaAccessor<structdim_, dim, spacedim>>(&tria, 0, pos);
    }
  } // namespace TriangulationImplementation
} // namespace internal
 
 
 
DEAL_II_NAMESPACE_CLOSE
 
#endif