Reference documentation for deal.II version GIT d77e5ebb0a 2023-01-27 22:35:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
tria_accessor.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2022 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_tria_accessor_h
17 #define dealii_tria_accessor_h
18 
19 
20 #include <deal.II/base/config.h>
21 
25 #include <deal.II/base/point.h>
26 
27 #include <deal.II/grid/cell_id.h>
31 
33 #include <boost/container/small_vector.hpp>
35 
36 #include <utility>
37 
38 
40 
41 // Forward declarations
42 #ifndef DOXYGEN
43 template <int dim, int spacedim>
44 class Triangulation;
45 template <typename Accessor>
46 class TriaRawIterator;
47 template <typename Accessor>
48 class TriaIterator;
49 template <typename Accessor>
50 class TriaActiveIterator;
51 
52 namespace parallel
53 {
54  template <int dim, int spacedim>
55  class TriangulationBase;
56 }
57 
58 template <int dim, int spacedim>
59 class DoFHandler;
60 template <int dim, int spacedim, bool lda>
61 class DoFCellAccessor;
62 
63 
64 template <int dim, int spacedim>
65 class Manifold;
66 
67 template <int dim, int spacedim>
68 class Mapping;
69 #endif
70 
71 namespace internal
72 {
73  namespace TriangulationImplementation
74  {
75  class TriaObjects;
76  struct Implementation;
77  struct ImplementationMixedMesh;
78  } // namespace TriangulationImplementation
79 
80  namespace TriaAccessorImplementation
81  {
82  struct Implementation;
83 
89  template <int structdim, int dim>
91  {
92  struct type
93  {
97  type() = default;
98 
102  type(const int level)
103  {
104  Assert(level == 0, ExcInternalError());
105  (void)level; // removes -Wunused-parameter warning in optimized mode
106  }
107 
111  operator int() const
112  {
113  return 0;
114  }
115 
116  void
117  operator++() const
118  {
119  Assert(false, ExcInternalError());
120  }
121 
122  void
123  operator--() const
124  {
125  Assert(false, ExcInternalError());
126  }
127  };
128  };
129 
130 
136  template <int dim>
137  struct PresentLevelType<dim, dim>
138  {
139  using type = int;
140  };
141  } // namespace TriaAccessorImplementation
142 } // namespace internal
143 template <int structdim, int dim, int spacedim>
144 class TriaAccessor;
145 template <int dim, int spacedim>
146 class TriaAccessor<0, dim, spacedim>;
147 template <int spacedim>
148 class TriaAccessor<0, 1, spacedim>;
149 
150 // note: the file tria_accessor.templates.h is included at the end of
151 // this file. this includes a lot of templates. originally, this was
152 // only done in debug mode, but led to cyclic reduction problems and
153 // so is now on by default.
154 
155 
160 {
165  "The operation you are attempting can only be performed for "
166  "(cell, face, or edge) iterators that point to valid "
167  "objects. These objects need not necessarily be active, "
168  "i.e., have no children, but they need to be part of a "
169  "triangulation. (The objects pointed to by an iterator "
170  "may -- after coarsening -- also be objects that used "
171  "to be part of a triangulation, but are now no longer "
172  "used. Their memory location may have been retained "
173  "for re-use upon the next mesh refinement, but is "
174  "currently unused.)");
185  "The operation you are attempting can only be performed for "
186  "(cell, face, or edge) iterators that point to 'active' "
187  "objects. 'Active' objects are those that do not have "
188  "children (in the case of cells), or that are part of "
189  "an active cell (in the case of faces or edges). However, "
190  "the object on which you are trying the current "
191  "operation is not 'active' in this sense.");
198  "The operation you are attempting can only be performed for "
199  "(cell, face, or edge) iterators that have children, "
200  "but the object on which you are trying the current "
201  "operation does not have any.");
209  "The operation you are attempting can only be performed for "
210  "(cell, face, or edge) iterators that have a parent object, "
211  "but the object on which you are trying the current "
212  "operation does not have one -- i.e., it is on the "
213  "coarsest level of the triangulation.");
218  int,
219  << "You can only set the child index if the cell does not "
220  << "currently have children registered; or you can clear it. "
221  << "The given index was " << arg1
222  << " (-1 means: clear children).");
226  template <typename AccessorType>
228  AccessorType,
229  << "You tried to dereference an iterator for which this "
230  << "is not possible. More information on this iterator: "
231  << "index=" << arg1.index() << ", state="
232  << (arg1.state() == IteratorState::valid ?
233  "valid" :
234  (arg1.state() == IteratorState::past_the_end ?
235  "past_the_end" :
236  "invalid")));
241  "Iterators can only be compared if they point to the same "
242  "triangulation, or if neither of them are associated "
243  "with a triangulation.");
244  // TODO: Write documentation!
249  // TODO: Write documentation!
269  // TODO: Write documentation!
275  int,
276  << "You can only set the child index of an even numbered child."
277  << "The number of the child given was " << arg1 << '.');
278 } // namespace TriaAccessorExceptions
279 
280 
306 template <int structdim, int dim, int spacedim = dim>
308 {
309 public:
315  static constexpr unsigned int space_dimension = spacedim;
316 
322  static constexpr unsigned int dimension = dim;
323 
329  static const unsigned int structure_dimension = structdim;
330 
340  void
341  operator=(const TriaAccessorBase *) = delete;
342 
343 protected:
349  using AccessorData = void;
350 
355  const int level = -1,
356  const int index = -1,
357  const AccessorData * = nullptr);
358 
363 
371  void
373 
379 
390  bool
391  operator<(const TriaAccessorBase &other) const;
392 
393 protected:
397  bool
398  operator==(const TriaAccessorBase &) const;
399 
403  bool
404  operator!=(const TriaAccessorBase &) const;
405 
419  void
421 
429  void
439  objects() const;
440 
441 public:
447  using LocalData = void *;
448 
472  int
473  level() const;
474 
501  int
502  index() const;
503 
509  state() const;
510 
517 
521 protected:
526  typename ::internal::TriaAccessorImplementation::
527  PresentLevelType<structdim, dim>::type present_level;
528 
534 
539 
540 private:
541  template <typename Accessor>
542  friend class TriaRawIterator;
543  template <typename Accessor>
544  friend class TriaIterator;
545  template <typename Accessor>
546  friend class TriaActiveIterator;
547 };
548 
549 
550 
571 template <int structdim, int dim, int spacedim = dim>
572 class InvalidAccessor : public TriaAccessorBase<structdim, dim, spacedim>
573 {
574 public:
578  using AccessorData =
580 
589  const int level = -1,
590  const int index = -1,
591  const AccessorData * local_data = nullptr);
592 
601 
606  template <typename OtherAccessor>
607  InvalidAccessor(const OtherAccessor &);
608 
612  void
614 
618  bool
619  operator==(const InvalidAccessor &) const;
620  bool
621  operator!=(const InvalidAccessor &) const;
622 
626  void
627  operator++() const;
628  void
629  operator--() const;
630 
635  bool
636  used() const;
637 
642  bool
643  has_children() const;
644 
649  manifold_id() const;
650 
654  unsigned int
655  user_index() const;
656 
660  void
661  set_user_index(const unsigned int p) const;
662 
666  void
668 
673  vertex(const unsigned int i) const;
674 
679  typename ::internal::TriangulationImplementation::
680  Iterators<dim, spacedim>::line_iterator
681  line(const unsigned int i) const;
682 
687  typename ::internal::TriangulationImplementation::
688  Iterators<dim, spacedim>::quad_iterator
689  quad(const unsigned int i) const;
690 };
691 
692 
693 
711 template <int structdim, int dim, int spacedim>
712 class TriaAccessor : public TriaAccessorBase<structdim, dim, spacedim>
713 {
714 public:
718  using AccessorData =
720 
725  const int level = -1,
726  const int index = -1,
727  const AccessorData * local_data = nullptr);
728 
733  TriaAccessor(const TriaAccessor &) = default;
734 
738  TriaAccessor(TriaAccessor &&) = default; // NOLINT
739 
752  template <int structdim2, int dim2, int spacedim2>
754 
759  template <int structdim2, int dim2, int spacedim2>
761 
771  TriaAccessor &
772  operator=(const TriaAccessor &) = delete;
773 
777  TriaAccessor &
778  operator=(TriaAccessor &&) = default; // NOLINT
779 
783  ~TriaAccessor() = default;
784 
791  bool
792  used() const;
793 
806  vertex_iterator(const unsigned int i) const;
807 
823  unsigned int
824  vertex_index(const unsigned int i) const;
825 
864  vertex(const unsigned int i) const;
865 
869  typename ::internal::TriangulationImplementation::
870  Iterators<dim, spacedim>::line_iterator
871  line(const unsigned int i) const;
872 
879  unsigned int
880  line_index(const unsigned int i) const;
881 
885  typename ::internal::TriangulationImplementation::
886  Iterators<dim, spacedim>::quad_iterator
887  quad(const unsigned int i) const;
888 
895  unsigned int
896  quad_index(const unsigned int i) const;
912  unsigned char
913  combined_face_orientation(const unsigned int face) const;
914 
926  bool
927  face_orientation(const unsigned int face) const;
928 
938  bool
939  face_flip(const unsigned int face) const;
940 
950  bool
951  face_rotation(const unsigned int face) const;
952 
965  bool
966  line_orientation(const unsigned int line) const;
981  bool
982  has_children() const;
983 
988  unsigned int
989  n_children() const;
990 
995  unsigned int
997 
1011  unsigned int
1013 
1027  unsigned int
1029 
1034  child(const unsigned int i) const;
1035 
1040  unsigned int
1043 
1053  isotropic_child(const unsigned int i) const;
1054 
1060 
1066  int
1067  child_index(const unsigned int i) const;
1068 
1074  int
1075  isotropic_child_index(const unsigned int i) const;
1098  boundary_id() const;
1099 
1129  void
1131 
1162  void
1164 
1172  bool
1173  at_boundary() const;
1174 
1184  const Manifold<dim, spacedim> &
1185  get_manifold() const;
1186 
1208  manifold_id() const;
1209 
1227  void
1229 
1243  void
1245 
1262  bool
1263  user_flag_set() const;
1264 
1270  void
1271  set_user_flag() const;
1272 
1278  void
1280 
1286  void
1288 
1294  void
1296 
1302  void
1304 
1316  void
1317  set_user_pointer(void *p) const;
1318 
1324  void
1326 
1342  void *
1343  user_pointer() const;
1344 
1366  void
1368 
1375  void
1377 
1387  void
1388  set_user_index(const unsigned int p) const;
1389 
1395  void
1397 
1409  unsigned int
1410  user_index() const;
1411 
1429  void
1430  recursively_set_user_index(const unsigned int p) const;
1431 
1440  void
1476  double
1477  diameter() const;
1478 
1505  std::pair<Point<spacedim>, double>
1507 
1517  bounding_box() const;
1518 
1528  double
1529  extent_in_direction(const unsigned int axis) const;
1530 
1534  double
1536 
1551  intermediate_point(const Point<structdim> &coordinates) const;
1552 
1577 
1613  center(const bool respect_manifold = false,
1614  const bool interpolate_from_surrounding = false) const;
1615 
1634  barycenter() const;
1635 
1661  double
1662  measure() const;
1663 
1678  bool
1681 
1687 
1691  unsigned int
1692  n_vertices() const;
1693 
1697  unsigned int
1698  n_lines() const;
1699 
1705  unsigned int
1706  n_faces() const;
1707 
1714 
1720  line_indices() const;
1721 
1729  face_indices() const;
1730 
1735 private:
1740  void
1742 
1750  void
1752  const std::initializer_list<int> &new_indices) const;
1753 
1757  void
1759  const std::initializer_list<unsigned int> &new_indices) const;
1760 
1768  void
1769  set_line_orientation(const unsigned int line, const bool orientation) const;
1770 
1778  void
1779  set_combined_face_orientation(const unsigned int face,
1780  const unsigned char combined_orientation) const;
1781 
1785  void
1786  set_used_flag() const;
1787 
1791  void
1793 
1802  void
1804 
1812  void
1814 
1821  void
1822  set_children(const unsigned int i, const int index) const;
1823 
1828  void
1830 
1831 private:
1832  template <int, int>
1833  friend class Triangulation;
1834 
1835  friend struct ::internal::TriangulationImplementation::Implementation;
1836  friend struct ::internal::TriangulationImplementation::
1837  ImplementationMixedMesh;
1838  friend struct ::internal::TriaAccessorImplementation::Implementation;
1839 };
1840 
1841 
1842 
1861 template <int dim, int spacedim>
1862 class TriaAccessor<0, dim, spacedim>
1863 {
1864 public:
1870  static constexpr unsigned int space_dimension = spacedim;
1871 
1877  static constexpr unsigned int dimension = dim;
1878 
1884  static const unsigned int structure_dimension = 0;
1885 
1889  using AccessorData = void;
1890 
1896  const unsigned int vertex_index);
1897 
1904  const int level = 0,
1905  const int index = 0,
1906  const AccessorData * = nullptr);
1907 
1911  template <int structdim2, int dim2, int spacedim2>
1913 
1917  template <int structdim2, int dim2, int spacedim2>
1919 
1924  state() const;
1925 
1930  static int
1932 
1937  int
1938  index() const;
1939 
1946 
1956  void
1958 
1962  void
1967  bool
1968  operator==(const TriaAccessor &) const;
1969 
1973  bool
1974  operator!=(const TriaAccessor &) const;
1975 
2003  unsigned int
2004  vertex_index(const unsigned int i = 0) const;
2005 
2011  Point<spacedim> &
2012  vertex(const unsigned int i = 0) const;
2013 
2018  typename ::internal::TriangulationImplementation::
2019  Iterators<dim, spacedim>::line_iterator static line(const unsigned int);
2020 
2024  static unsigned int
2025  line_index(const unsigned int i);
2026 
2030  static typename ::internal::TriangulationImplementation::
2031  Iterators<dim, spacedim>::quad_iterator
2032  quad(const unsigned int i);
2033 
2037  static unsigned int
2038  quad_index(const unsigned int i);
2039 
2055  double
2056  diameter() const;
2057 
2065  double
2066  extent_in_direction(const unsigned int axis) const;
2067 
2076  center(const bool respect_manifold = false,
2077  const bool interpolate_from_surrounding = false) const;
2078 
2087  double
2088  measure() const;
2103  static unsigned char
2104  combined_face_orientation(const unsigned int face);
2105 
2109  static bool
2110  face_orientation(const unsigned int face);
2111 
2115  static bool
2116  face_flip(const unsigned int face);
2117 
2121  static bool
2122  face_rotation(const unsigned int face);
2123 
2127  static bool
2128  line_orientation(const unsigned int line);
2129 
2144  static bool
2146 
2151  static unsigned int
2153 
2158  static unsigned int
2160 
2165  static unsigned int
2167 
2168 
2172  static unsigned int
2174 
2178  static unsigned int
2180 
2185  child(const unsigned int);
2186 
2191  isotropic_child(const unsigned int);
2192 
2196  static RefinementCase<0>
2198 
2202  static int
2203  child_index(const unsigned int i);
2204 
2208  static int
2209  isotropic_child_index(const unsigned int i);
2217  bool
2218  used() const;
2219 
2220 protected:
2228  void
2230 
2239  bool
2240  operator<(const TriaAccessor &other) const;
2241 
2246 
2250  unsigned int global_vertex_index;
2251 
2252 private:
2253  template <typename Accessor>
2254  friend class TriaRawIterator;
2255  template <typename Accessor>
2256  friend class TriaIterator;
2257  template <typename Accessor>
2258  friend class TriaActiveIterator;
2259 };
2260 
2261 
2262 
2279 template <int spacedim>
2280 class TriaAccessor<0, 1, spacedim>
2281 {
2282 public:
2288  static constexpr unsigned int space_dimension = spacedim;
2289 
2295  static constexpr unsigned int dimension = 1;
2296 
2302  static const unsigned int structure_dimension = 0;
2303 
2307  using AccessorData = void;
2308 
2314  {
2326  right_vertex
2327  };
2328 
2341  const VertexKind vertex_kind,
2342  const unsigned int vertex_index);
2343 
2350  const int = 0,
2351  const int = 0,
2352  const AccessorData * = nullptr);
2353 
2357  template <int structdim2, int dim2, int spacedim2>
2359 
2363  template <int structdim2, int dim2, int spacedim2>
2365 
2370  void
2372 
2378  void
2380 
2388 
2393  static int
2395 
2400  int
2401  index() const;
2402 
2409 
2420  void
2421  operator++() const;
2422 
2427  void
2428  operator--() const;
2432  bool
2433  operator==(const TriaAccessor &) const;
2434 
2438  bool
2439  operator!=(const TriaAccessor &) const;
2440 
2449  bool
2450  operator<(const TriaAccessor &other) const;
2451 
2478  unsigned int
2479  vertex_index(const unsigned int i = 0) const;
2480 
2486  Point<spacedim> &
2487  vertex(const unsigned int i = 0) const;
2488 
2494  center() const;
2495 
2500  typename ::internal::TriangulationImplementation::
2501  Iterators<1, spacedim>::line_iterator static line(const unsigned int);
2502 
2509  static unsigned int
2510  line_index(const unsigned int i);
2511 
2515  static typename ::internal::TriangulationImplementation::
2516  Iterators<1, spacedim>::quad_iterator
2517  quad(const unsigned int i);
2518 
2525  static unsigned int
2526  quad_index(const unsigned int i);
2527 
2537  bool
2538  at_boundary() const;
2539 
2555  boundary_id() const;
2556 
2560  const Manifold<1, spacedim> &
2561  get_manifold() const;
2562 
2570  manifold_id() const;
2571 
2572 
2584  bool
2585  user_flag_set() const;
2586 
2592  void
2593  set_user_flag() const;
2594 
2600  void
2601  clear_user_flag() const;
2602 
2608  void
2609  recursively_set_user_flag() const;
2610 
2616  void
2618 
2624  void
2625  clear_user_data() const;
2626 
2638  void
2639  set_user_pointer(void *p) const;
2640 
2646  void
2647  clear_user_pointer() const;
2648 
2664  void *
2665  user_pointer() const;
2666 
2688  void
2689  recursively_set_user_pointer(void *p) const;
2690 
2697  void
2699 
2709  void
2710  set_user_index(const unsigned int p) const;
2711 
2717  void
2718  clear_user_index() const;
2719 
2731  unsigned int
2732  user_index() const;
2733 
2751  void
2752  recursively_set_user_index(const unsigned int p) const;
2753 
2762  void
2778  static unsigned char
2779  combined_face_orientation(const unsigned int face);
2780 
2784  static bool
2785  face_orientation(const unsigned int face);
2786 
2790  static bool
2791  face_flip(const unsigned int face);
2792 
2796  static bool
2797  face_rotation(const unsigned int face);
2798 
2802  static bool
2803  line_orientation(const unsigned int line);
2804 
2819  static bool
2821 
2826  static unsigned int
2828 
2833  static unsigned int
2835 
2840  static unsigned int
2842 
2843 
2847  static unsigned int
2849 
2853  static unsigned int
2855 
2860  child(const unsigned int);
2861 
2866  isotropic_child(const unsigned int);
2867 
2871  static RefinementCase<0>
2873 
2877  static int
2878  child_index(const unsigned int i);
2879 
2883  static int
2884  isotropic_child_index(const unsigned int i);
2917  void
2919 
2926  void
2928 
2940  void
2942 
2954  void
2963  bool
2964  used() const;
2965 
2971 
2975  unsigned int
2976  n_vertices() const;
2977 
2981  unsigned int
2982  n_lines() const;
2983 
2990 
2996  line_indices() const;
2997 
2998 protected:
3003 
3009 
3013  unsigned int global_vertex_index;
3014 };
3015 
3016 
3017 
3033 template <int dim, int spacedim = dim>
3034 class CellAccessor : public TriaAccessor<dim, dim, spacedim>
3035 {
3036 public:
3041 
3046 
3058  const int level = -1,
3059  const int index = -1,
3060  const AccessorData * local_data = nullptr);
3061 
3066 
3079  template <int structdim2, int dim2, int spacedim2>
3081 
3086  template <int structdim2, int dim2, int spacedim2>
3088 
3093 
3097  // NOLINTNEXTLINE OSX does not compile with noexcept
3099 
3103  ~CellAccessor() = default;
3104 
3116 
3120  // NOLINTNEXTLINE OSX does not compile with noexcept
3122  operator=(CellAccessor<dim, spacedim> &&) = default; // NOLINT
3123 
3147  as_dof_handler_iterator(const DoFHandler<dim, spacedim> &dof_handler) const;
3148 
3165  child(const unsigned int i) const;
3166 
3170  boost::container::small_vector<TriaIterator<CellAccessor<dim, spacedim>>,
3173 
3177  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>
3178  face(const unsigned int i) const;
3179 
3184  unsigned int
3187 
3191  boost::container::small_vector<
3192  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>,
3195 
3205  unsigned int
3206  face_index(const unsigned int i) const;
3207 
3256  neighbor_child_on_subface(const unsigned int face_no,
3257  const unsigned int subface_no) const;
3258 
3284  neighbor(const unsigned int face_no) const;
3285 
3293  int
3294  neighbor_index(const unsigned int face_no) const;
3295 
3303  int
3304  neighbor_level(const unsigned int face_no) const;
3305 
3317  unsigned int
3318  neighbor_of_neighbor(const unsigned int face_no) const;
3319 
3330  bool
3331  neighbor_is_coarser(const unsigned int face_no) const;
3332 
3347  std::pair<unsigned int, unsigned int>
3348  neighbor_of_coarser_neighbor(const unsigned int neighbor) const;
3349 
3356  unsigned int
3357  neighbor_face_no(const unsigned int neighbor) const;
3358 
3362  static bool
3364 
3378  bool
3379  has_periodic_neighbor(const unsigned int i) const;
3380 
3398  periodic_neighbor(const unsigned int i) const;
3399 
3408  neighbor_or_periodic_neighbor(const unsigned int i) const;
3409 
3425  periodic_neighbor_child_on_subface(const unsigned int face_no,
3426  const unsigned int subface_no) const;
3427 
3438  std::pair<unsigned int, unsigned int>
3439  periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const;
3440 
3446  int
3447  periodic_neighbor_index(const unsigned int i) const;
3448 
3454  int
3455  periodic_neighbor_level(const unsigned int i) const;
3456 
3471  unsigned int
3472  periodic_neighbor_of_periodic_neighbor(const unsigned int i) const;
3473 
3479  unsigned int
3480  periodic_neighbor_face_no(const unsigned int i) const;
3481 
3488  bool
3489  periodic_neighbor_is_coarser(const unsigned int i) const;
3490 
3507  bool
3508  at_boundary(const unsigned int i) const;
3509 
3518  bool
3519  at_boundary() const;
3520 
3528  bool
3529  has_boundary_lines() const;
3557 
3575  void
3578 
3582  void
3584 
3592  bool
3594  const unsigned int face_no,
3595  const RefinementCase<dim - 1> &face_refinement_case =
3597 
3603  bool
3604  flag_for_line_refinement(const unsigned int line_no) const;
3605 
3615  subface_case(const unsigned int face_no) const;
3616 
3620  bool
3622 
3627  void
3629 
3633  void
3658  material_id() const;
3659 
3671  void
3672  set_material_id(const types::material_id new_material_id) const;
3673 
3682  void
3683  recursively_set_material_id(const types::material_id new_material_id) const;
3710  subdomain_id() const;
3711 
3727  void
3728  set_subdomain_id(const types::subdomain_id new_subdomain_id) const;
3729 
3740 
3745  void
3747  const types::subdomain_id new_level_subdomain_id) const;
3748 
3749 
3765  void
3767  const types::subdomain_id new_subdomain_id) const;
3791 
3801 
3815  bool
3816  direction_flag() const;
3817 
3843  unsigned int
3845 
3853  int
3854  parent_index() const;
3855 
3862  parent() const;
3863 
3883  bool
3884  is_active() const;
3885 
3905  bool
3907 
3912  bool
3914 
3948  bool
3949  is_ghost() const;
3950 
3956  bool
3958 
3985  bool
3986  is_artificial() const;
3987 
3994  bool
3996 
4010  bool
4012 
4021  void
4022  set_neighbor(const unsigned int i,
4023  const TriaIterator<CellAccessor<dim, spacedim>> &pointer) const;
4024 
4038  CellId
4039  id() const;
4040 
4042 
4046  double
4047  diameter(const Mapping<dim, spacedim> &mapping) const;
4048 
4066 
4067 protected:
4083  unsigned int
4084  neighbor_of_neighbor_internal(const unsigned int neighbor) const;
4085 
4091  template <int dim_, int spacedim_>
4092  bool
4093  point_inside_codim(const Point<spacedim_> &p) const;
4094 
4095 
4096 
4097 private:
4102  void
4103  set_active_cell_index(const unsigned int active_cell_index) const;
4104 
4108  void
4110 
4114  void
4116 
4120  void
4121  set_parent(const unsigned int parent_index);
4122 
4129  void
4130  set_direction_flag(const bool new_direction_flag) const;
4131 
4132  template <int, int>
4133  friend class Triangulation;
4134 
4135  template <int, int>
4137 
4138  friend struct ::internal::TriangulationImplementation::Implementation;
4139  friend struct ::internal::TriangulationImplementation::
4140  ImplementationMixedMesh;
4141 };
4142 
4143 
4144 
4145 /* ----- declaration of explicit specializations and general templates ----- */
4146 
4147 
4148 template <int structdim, int dim, int spacedim>
4149 template <typename OtherAccessor>
4151  const OtherAccessor &)
4152 {
4153  Assert(false,
4154  ExcMessage("You are attempting an illegal conversion between "
4155  "iterator/accessor types. The constructor you call "
4156  "only exists to make certain template constructs "
4157  "easier to write as dimension independent code but "
4158  "the conversion is not valid in the current context."));
4159 }
4160 
4161 
4162 
4163 template <int structdim, int dim, int spacedim>
4164 template <int structdim2, int dim2, int spacedim2>
4167 {
4168  Assert(false,
4169  ExcMessage("You are attempting an illegal conversion between "
4170  "iterator/accessor types. The constructor you call "
4171  "only exists to make certain template constructs "
4172  "easier to write as dimension independent code but "
4173  "the conversion is not valid in the current context."));
4174 }
4175 
4176 
4177 
4178 template <int dim, int spacedim>
4179 template <int structdim2, int dim2, int spacedim2>
4182 {
4183  Assert(false,
4184  ExcMessage("You are attempting an illegal conversion between "
4185  "iterator/accessor types. The constructor you call "
4186  "only exists to make certain template constructs "
4187  "easier to write as dimension independent code but "
4188  "the conversion is not valid in the current context."));
4189 }
4190 
4191 
4192 
4193 template <int structdim, int dim, int spacedim>
4194 template <int structdim2, int dim2, int spacedim2>
4197 {
4198  Assert(false,
4199  ExcMessage("You are attempting an illegal conversion between "
4200  "iterator/accessor types. The constructor you call "
4201  "only exists to make certain template constructs "
4202  "easier to write as dimension independent code but "
4203  "the conversion is not valid in the current context."));
4204 }
4205 
4206 
4207 
4208 template <int dim, int spacedim>
4209 template <int structdim2, int dim2, int spacedim2>
4212 {
4213  Assert(false,
4214  ExcMessage("You are attempting an illegal conversion between "
4215  "iterator/accessor types. The constructor you call "
4216  "only exists to make certain template constructs "
4217  "easier to write as dimension independent code but "
4218  "the conversion is not valid in the current context."));
4219 }
4220 
4221 
4222 #ifndef DOXYGEN
4223 
4224 template <>
4225 bool
4227 template <>
4228 bool
4230 template <>
4231 bool
4233 template <>
4234 bool
4236 template <>
4237 bool
4239 template <>
4240 bool
4242 // -------------------------------------------------------------------
4243 
4244 template <>
4245 void
4247 
4248 #endif // DOXYGEN
4249 
4251 
4252 // include more templates in debug and optimized mode
4253 #include "tria_accessor.templates.h"
4254 
4255 #endif
void recursively_set_subdomain_id(const types::subdomain_id new_subdomain_id) const
CellAccessor(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
TriaIterator< CellAccessor< dim, spacedim > > parent() const
unsigned int neighbor_face_no(const unsigned int neighbor) const
unsigned int face_iterator_to_index(const TriaIterator< TriaAccessor< dim - 1, dim, spacedim >> &face) const
void set_active_cell_index(const unsigned int active_cell_index) const
unsigned int neighbor_of_neighbor_internal(const unsigned int neighbor) const
TriaIterator< CellAccessor< dim, spacedim > > periodic_neighbor(const unsigned int i) const
types::global_cell_index global_active_cell_index() const
bool is_artificial_on_level() const
void set_direction_flag(const bool new_direction_flag) const
void recursively_set_material_id(const types::material_id new_material_id) const
void set_level_subdomain_id(const types::subdomain_id new_level_subdomain_id) const
types::subdomain_id level_subdomain_id() const
boost::container::small_vector< TriaIterator< CellAccessor< dim, spacedim > >, GeometryInfo< dim >::max_children_per_cell > child_iterators() const
bool flag_for_face_refinement(const unsigned int face_no, const RefinementCase< dim - 1 > &face_refinement_case=RefinementCase< dim - 1 >::isotropic_refinement) const
unsigned int face_index(const unsigned int i) const
double diameter(const Mapping< dim, spacedim > &mapping) const
TriaActiveIterator< DoFCellAccessor< dim, spacedim, false > > as_dof_handler_iterator(const DoFHandler< dim, spacedim > &dof_handler) const
::internal::SubfaceCase< dim > subface_case(const unsigned int face_no) const
bool is_active() const
bool is_ghost() const
TriaIterator< CellAccessor< dim, spacedim > > neighbor_child_on_subface(const unsigned int face_no, const unsigned int subface_no) const
void set_subdomain_id(const types::subdomain_id new_subdomain_id) const
CellAccessor< dim, spacedim > & operator=(CellAccessor< dim, spacedim > &&)=default
bool neighbor_is_coarser(const unsigned int face_no) const
void set_global_level_cell_index(const types::global_cell_index index) const
bool has_periodic_neighbor(const unsigned int i) const
int periodic_neighbor_level(const unsigned int i) const
std::pair< unsigned int, unsigned int > neighbor_of_coarser_neighbor(const unsigned int neighbor) const
~CellAccessor()=default
CellAccessor(const CellAccessor< dim, spacedim > &)=default
void set_coarsen_flag() const
TriaIterator< TriaAccessor< dim - 1, dim, spacedim > > face(const unsigned int i) const
unsigned int neighbor_of_neighbor(const unsigned int face_no) const
void set_material_id(const types::material_id new_material_id) const
bool is_locally_owned() const
TriaIterator< CellAccessor< dim, spacedim > > neighbor(const unsigned int face_no) const
void set_refine_flag(const RefinementCase< dim > ref_case=RefinementCase< dim >::isotropic_refinement) const
RefinementCase< dim > refine_flag_set() const
CellAccessor(CellAccessor< dim, spacedim > &&)=default
bool point_inside_codim(const Point< spacedim_ > &p) const
typename TriaAccessor< dim, dim, spacedim >::AccessorData AccessorData
bool is_locally_owned_on_level() const
bool has_boundary_lines() const
TriaIterator< CellAccessor< dim, spacedim > > periodic_neighbor_child_on_subface(const unsigned int face_no, const unsigned int subface_no) const
int neighbor_level(const unsigned int face_no) const
int periodic_neighbor_index(const unsigned int i) const
bool periodic_neighbor_is_coarser(const unsigned int i) const
void set_global_active_cell_index(const types::global_cell_index index) const
void clear_coarsen_flag() const
TriaIterator< CellAccessor< dim, spacedim > > child(const unsigned int i) const
void set_neighbor(const unsigned int i, const TriaIterator< CellAccessor< dim, spacedim >> &pointer) const
void set_parent(const unsigned int parent_index)
std::pair< unsigned int, unsigned int > periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const
CellAccessor(const TriaAccessor< dim, dim, spacedim > &cell_accessor)
bool at_boundary() const
unsigned int active_cell_index() const
static bool is_level_cell()
void clear_refine_flag() const
int neighbor_index(const unsigned int face_no) const
bool point_inside(const Point< spacedim > &p) const
types::subdomain_id subdomain_id() const
bool direction_flag() const
types::material_id material_id() const
bool coarsen_flag_set() const
types::global_cell_index global_level_cell_index() const
bool flag_for_line_refinement(const unsigned int line_no) const
boost::container::small_vector< TriaIterator< TriaAccessor< dim - 1, dim, spacedim > >, GeometryInfo< dim >::faces_per_cell > face_iterators() const
CellId id() const
bool is_artificial() const
bool is_ghost_on_level() const
TriaIterator< CellAccessor< dim, spacedim > > neighbor_or_periodic_neighbor(const unsigned int i) const
int parent_index() const
unsigned int periodic_neighbor_of_periodic_neighbor(const unsigned int i) const
unsigned int periodic_neighbor_face_no(const unsigned int i) const
CellAccessor< dim, spacedim > & operator=(const CellAccessor< dim, spacedim > &)=delete
Definition: cell_id.h:71
typename ::internal::TriangulationImplementation::Iterators< dim, spacedim >::quad_iterator quad(const unsigned int i) const
void set_manifold_id(const types::manifold_id) const
InvalidAccessor(const InvalidAccessor &)
typename ::internal::TriangulationImplementation::Iterators< dim, spacedim >::line_iterator line(const unsigned int i) const
InvalidAccessor(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
bool operator==(const InvalidAccessor &) const
void operator++() const
bool operator!=(const InvalidAccessor &) const
bool has_children() const
unsigned int user_index() const
Point< spacedim > & vertex(const unsigned int i) const
void copy_from(const InvalidAccessor &)
void set_user_index(const unsigned int p) const
types::manifold_id manifold_id() const
void operator--() const
bool used() const
InvalidAccessor(const OtherAccessor &)
typename TriaAccessorBase< structdim, dim, spacedim >::AccessorData AccessorData
Abstract base class for mapping classes.
Definition: mapping.h:312
bool operator!=(const TriaAccessorBase &) const
static constexpr unsigned int space_dimension
static constexpr unsigned int dimension
TriaAccessorBase(const TriaAccessorBase &)
void operator=(const TriaAccessorBase *)=delete
static const unsigned int structure_dimension
typename ::internal::TriaAccessorImplementation::PresentLevelType< structdim, dim >::type present_level
TriaAccessorBase(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *=nullptr)
void copy_from(const TriaAccessorBase &)
IteratorState::IteratorStates state() const
int index() const
bool operator<(const TriaAccessorBase &other) const
TriaAccessorBase & operator=(const TriaAccessorBase &)
::internal::TriangulationImplementation::TriaObjects & objects() const
int level() const
const Triangulation< dim, spacedim > * tria
const Triangulation< dim, spacedim > & get_triangulation() const
bool operator==(const TriaAccessorBase &) const
static unsigned int line_index(const unsigned int i)
static bool face_orientation(const unsigned int face)
Always return false.
static unsigned int number_of_children()
static RefinementCase< 0 > refinement_case()
static unsigned int max_refinement_depth()
Point< spacedim > center() const
static unsigned int quad_index(const unsigned int i)
TriaAccessor(const InvalidAccessor< structdim2, dim2, spacedim2 > &)
static bool face_flip(const unsigned int face)
Always return false.
TriaAccessor(const Triangulation< 1, spacedim > *tria, const VertexKind vertex_kind, const unsigned int vertex_index)
static TriaIterator< TriaAccessor< 0, 1, spacedim > > isotropic_child(const unsigned int)
Return an invalid object.
std_cxx20::ranges::iota_view< unsigned int, unsigned int > line_indices() const
unsigned int n_lines() const
bool operator!=(const TriaAccessor &) const
TriaAccessor(const Triangulation< 1, spacedim > *tria=nullptr, const int=0, const int=0, const AccessorData *=nullptr)
static unsigned int n_active_descendants()
static int isotropic_child_index(const unsigned int i)
Returns -1.
types::boundary_id boundary_id() const
static typename ::internal::TriangulationImplementation::Iterators< 1, spacedim >::line_iterator line(const unsigned int)
Point< spacedim > & vertex(const unsigned int i=0) const
static bool line_orientation(const unsigned int line)
Always return false.
static typename ::internal::TriangulationImplementation::Iterators< 1, spacedim >::quad_iterator quad(const unsigned int i)
ReferenceCell reference_cell() const
static unsigned int child_iterator_to_index(const TriaIterator< TriaAccessor< 0, 1, spacedim >> &)
Return an invalid unsigned integer.
unsigned int vertex_index(const unsigned int i=0) const
const Triangulation< 1, spacedim > & get_triangulation() const
void copy_from(const TriaAccessor &)
void set_all_boundary_ids(const types::boundary_id) const
static TriaIterator< TriaAccessor< 0, 1, spacedim > > child(const unsigned int)
Return an invalid object.
static IteratorState::IteratorStates state()
static bool face_rotation(const unsigned int face)
Always return false.
static int child_index(const unsigned int i)
Returns -1.
const Manifold< 1, spacedim > & get_manifold() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices() const
types::manifold_id manifold_id() const
static unsigned int n_children()
TriaAccessor(const TriaAccessor< structdim2, dim2, spacedim2 > &)
static unsigned char combined_face_orientation(const unsigned int face)
Always return 0.
void set_boundary_id(const types::boundary_id) const
unsigned int n_vertices() const
void set_manifold_id(const types::manifold_id)
void copy_from(const TriaAccessorBase< 0, 1, spacedim > &)
bool operator==(const TriaAccessor &) const
const Triangulation< 1, spacedim > * tria
double extent_in_direction(const unsigned int axis) const
TriaAccessor(const Triangulation< dim, spacedim > *tria, const unsigned int vertex_index)
static RefinementCase< 0 > refinement_case()
bool operator!=(const TriaAccessor &) const
static typename ::internal::TriangulationImplementation::Iterators< dim, spacedim >::line_iterator line(const unsigned int)
void copy_from(const TriaAccessor &)
static typename ::internal::TriangulationImplementation::Iterators< dim, spacedim >::quad_iterator quad(const unsigned int i)
TriaAccessor(const InvalidAccessor< structdim2, dim2, spacedim2 > &)
static int child_index(const unsigned int i)
Returns -1.
bool operator==(const TriaAccessor &) const
static unsigned int line_index(const unsigned int i)
unsigned int vertex_index(const unsigned int i=0) const
const Triangulation< dim, spacedim > * tria
static int isotropic_child_index(const unsigned int i)
Returns -1.
Point< spacedim > center(const bool respect_manifold=false, const bool interpolate_from_surrounding=false) const
TriaAccessor(const Triangulation< dim, spacedim > *tria=nullptr, const int level=0, const int index=0, const AccessorData *=nullptr)
static unsigned char combined_face_orientation(const unsigned int face)
Always return 0.
static unsigned int n_children()
static bool face_flip(const unsigned int face)
Always return false.
static TriaIterator< TriaAccessor< 0, dim, spacedim > > isotropic_child(const unsigned int)
Return an invalid object.
static bool line_orientation(const unsigned int line)
Always return false.
static unsigned int child_iterator_to_index(const TriaIterator< TriaAccessor< 0, dim, spacedim >> &)
Return an invalid unsigned integer.
static unsigned int max_refinement_depth()
static unsigned int quad_index(const unsigned int i)
static bool face_rotation(const unsigned int face)
Always return false.
const Triangulation< dim, spacedim > & get_triangulation() const
static bool face_orientation(const unsigned int face)
Always return false.
static unsigned int n_active_descendants()
static TriaIterator< TriaAccessor< 0, dim, spacedim > > child(const unsigned int)
Return an invalid object.
TriaAccessor(const TriaAccessor< structdim2, dim2, spacedim2 > &)
IteratorState::IteratorStates state() const
Point< spacedim > & vertex(const unsigned int i=0) const
static unsigned int number_of_children()
void clear_children() const
typename ::internal::TriangulationImplementation::Iterators< dim, spacedim >::quad_iterator quad(const unsigned int i) const
bool line_orientation(const unsigned int line) const
void set_boundary_id_internal(const types::boundary_id id) const
void set_user_index(const unsigned int p) const
TriaAccessor(TriaAccessor &&)=default
void clear_user_pointer() const
unsigned int child_iterator_to_index(const TriaIterator< TriaAccessor< structdim, dim, spacedim >> &child) const
unsigned int n_active_descendants() const
unsigned char combined_face_orientation(const unsigned int face) const
bool face_rotation(const unsigned int face) const
void recursively_set_user_index(const unsigned int p) const
void clear_user_data() const
TriaAccessor(const TriaAccessor< structdim2, dim2, spacedim2 > &)
Point< spacedim > & vertex(const unsigned int i) const
Point< structdim > real_to_unit_cell_affine_approximation(const Point< spacedim > &point) const
void recursively_clear_user_index() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > line_indices() const
unsigned int line_index(const unsigned int i) const
bool face_orientation(const unsigned int face) const
void recursively_set_user_pointer(void *p) const
unsigned int n_lines() const
double extent_in_direction(const unsigned int axis) const
Point< spacedim > intermediate_point(const Point< structdim > &coordinates) const
unsigned int n_vertices() const
bool has_children() const
void recursively_clear_user_flag() const
typename ::internal::TriangulationImplementation::Iterators< dim, spacedim >::line_iterator line(const unsigned int i) const
typename TriaAccessorBase< structdim, dim, spacedim >::AccessorData AccessorData
Point< spacedim > barycenter() const
BoundingBox< spacedim > bounding_box() const
void clear_user_flag() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > face_indices() const
TriaAccessor & operator=(TriaAccessor &&)=default
TriaAccessor & operator=(const TriaAccessor &)=delete
unsigned int n_children() const
void recursively_set_user_flag() const
void set_boundary_id(const types::boundary_id) const
bool user_flag_set() const
void set_used_flag() const
types::manifold_id manifold_id() const
void set_user_flag() const
TriaAccessor(const TriaAccessor &)=default
TriaIterator< TriaAccessor< 0, dim, spacedim > > vertex_iterator(const unsigned int i) const
void clear_used_flag() const
std::pair< Point< spacedim >, double > enclosing_ball() const
unsigned int vertex_index(const unsigned int i) const
int isotropic_child_index(const unsigned int i) const
Point< spacedim > center(const bool respect_manifold=false, const bool interpolate_from_surrounding=false) const
~TriaAccessor()=default
void set_refinement_case(const RefinementCase< structdim > &ref_case) const
void clear_user_index() const
double minimum_vertex_distance() const
double measure() const
void set_all_boundary_ids(const types::boundary_id) const
TriaIterator< TriaAccessor< structdim, dim, spacedim > > isotropic_child(const unsigned int i) const
void set_bounding_object_indices(const std::initializer_list< int > &new_indices) const
unsigned int number_of_children() const
unsigned int max_refinement_depth() const
void set_line_orientation(const unsigned int line, const bool orientation) const
bool is_translation_of(const TriaIterator< TriaAccessor< structdim, dim, spacedim >> &o) const
unsigned int quad_index(const unsigned int i) const
unsigned int user_index() const
int child_index(const unsigned int i) const
void set_user_pointer(void *p) const
void recursively_clear_user_pointer() const
ReferenceCell reference_cell() const
void * user_pointer() const
TriaIterator< TriaAccessor< structdim, dim, spacedim > > child(const unsigned int i) const
bool face_flip(const unsigned int face) const
TriaAccessor(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices() const
const Manifold< dim, spacedim > & get_manifold() const
RefinementCase< structdim > refinement_case() const
TriaAccessor(const InvalidAccessor< structdim2, dim2, spacedim2 > &)
bool used() const
void set_children(const unsigned int i, const int index) const
void clear_refinement_case() const
double diameter() const
bool at_boundary() const
types::boundary_id boundary_id() const
unsigned int n_faces() const
void set_combined_face_orientation(const unsigned int face, const unsigned char combined_orientation) const
#define DEAL_II_DEPRECATED
Definition: config.h:162
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:461
#define DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
Definition: config.h:475
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:462
#define DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
Definition: config.h:516
unsigned int level
Definition: grid_out.cc:4608
#define DeclException0(Exception0)
Definition: exceptions.h:464
static ::ExceptionBase & ExcRefineCellNotActive()
static ::ExceptionBase & ExcCantSetChildren(int arg1)
static ::ExceptionBase & ExcCellFlaggedForCoarsening()
static ::ExceptionBase & ExcCellNotActive()
static ::ExceptionBase & ExcCellHasNoParent()
static ::ExceptionBase & ExcNeighborIsNotCoarser()
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcCantCompareIterators()
static ::ExceptionBase & ExcCellNotUsed()
#define Assert(cond, exc)
Definition: exceptions.h:1509
static ::ExceptionBase & ExcCellFlaggedForRefinement()
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:487
static ::ExceptionBase & ExcDereferenceInvalidObject(AccessorType arg1)
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:509
static ::ExceptionBase & ExcSetOnlyEvenChildren(int arg1)
static ::ExceptionBase & ExcCellHasNoChildren()
static ::ExceptionBase & ExcMessage(std::string arg1)
static ::ExceptionBase & ExcFacesHaveNoLevel()
static ::ExceptionBase & ExcNoPeriodicNeighbor()
static ::ExceptionBase & ExcNeighborIsCoarser()
void set_all_manifold_ids(const types::manifold_id) const
void set_all_manifold_ids(const types::manifold_id)
void set_manifold_id(const types::manifold_id) const
@ past_the_end
Iterator reached end of container.
@ valid
Iterator points to a valid object.
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
Definition: utilities.cc:189
boost::integer_range< IncrementableType > iota_view
Definition: iota_view.h:46
unsigned int manifold_id
Definition: types.h:146
unsigned int subdomain_id
Definition: types.h:43
unsigned int global_cell_index
Definition: types.h:110
unsigned int material_id
Definition: types.h:157
unsigned int boundary_id
Definition: types.h:134