Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
tria_accessor.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2019 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 
22 #include <deal.II/base/bounding_box.h>
23 #include <deal.II/base/exceptions.h>
24 #include <deal.II/base/geometry_info.h>
25 #include <deal.II/base/point.h>
26 
27 #include <deal.II/grid/cell_id.h>
28 #include <deal.II/grid/tria_iterator_base.h>
29 #include <deal.II/grid/tria_iterator_selector.h>
30 
31 #include <utility>
32 
33 
34 DEAL_II_NAMESPACE_OPEN
35 
36 template <int dim, int spacedim>
37 class Triangulation;
38 template <typename Accessor>
39 class TriaRawIterator;
40 template <typename Accessor>
41 class TriaIterator;
42 template <typename Accessor>
43 class TriaActiveIterator;
44 
45 template <int dim, int spacedim>
46 class Manifold;
47 
48 
49 namespace internal
50 {
51  namespace TriangulationImplementation
52  {
53  template <int dim>
54  class TriaObject;
55  template <typename G>
56  class TriaObjects;
57  struct Implementation;
58  } // namespace TriangulationImplementation
59 
60  namespace TriaAccessorImplementation
61  {
62  struct Implementation;
63 
69  template <int structdim, int dim>
71  {
72  struct type
73  {
77  type() = default;
78 
82  type(const int level)
83  {
84  Assert(level == 0, ExcInternalError());
85  (void)level; // removes -Wunused-parameter warning in optimized mode
86  }
87 
91  operator int() const
92  {
93  return 0;
94  }
95 
96  void
97  operator++() const
98  {
99  Assert(false, ExcInternalError());
100  }
101 
102  void
103  operator--() const
104  {
105  Assert(false, ExcInternalError());
106  }
107  };
108  };
109 
110 
116  template <int dim>
117  struct PresentLevelType<dim, dim>
118  {
119  using type = int;
120  };
121 
122  } // namespace TriaAccessorImplementation
123 } // namespace internal
124 template <int structdim, int dim, int spacedim>
125 class TriaAccessor;
126 template <int dim, int spacedim>
127 class TriaAccessor<0, dim, spacedim>;
128 template <int spacedim>
129 class TriaAccessor<0, 1, spacedim>;
130 
131 // note: the file tria_accessor.templates.h is included at the end of
132 // this file. this includes a lot of templates. originally, this was
133 // only done in debug mode, but led to cyclic reduction problems and
134 // so is now on by default.
135 
136 
141 {
146  "The operation you are attempting can only be performed for "
147  "(cell, face, or edge) iterators that point to valid "
148  "objects. These objects need not necessarily be active, "
149  "i.e., have no children, but they need to be part of a "
150  "triangulation. (The objects pointed to by an iterator "
151  "may -- after coarsening -- also be objects that used "
152  "to be part of a triangulation, but are now no longer "
153  "used. Their memory location may have been retained "
154  "for re-use upon the next mesh refinement, but is "
155  "currently unused.)");
166  "The operation you are attempting can only be performed for "
167  "(cell, face, or edge) iterators that point to 'active' "
168  "objects. 'Active' objects are those that do not have "
169  "children (in the case of cells), or that are part of "
170  "an active cell (in the case of faces or edges). However, "
171  "the object on which you are trying the current "
172  "operation is not 'active' in this sense.");
179  "The operation you are attempting can only be performed for "
180  "(cell, face, or edge) iterators that have children, "
181  "but the object on which you are trying the current "
182  "operation does not have any.");
190  "The operation you are attempting can only be performed for "
191  "(cell, face, or edge) iterators that have a parent object, "
192  "but the object on which you are trying the current "
193  "operation does not have one -- i.e., it is on the "
194  "coarsest level of the triangulation.");
199  int,
200  << "You can only set the child index if the cell does not "
201  << "currently have children registered; or you can clear it. "
202  << "The given index was " << arg1
203  << " (-1 means: clear children).");
207  template <typename AccessorType>
209  AccessorType,
210  << "You tried to dereference an iterator for which this "
211  << "is not possible. More information on this iterator: "
212  << "index=" << arg1.index() << ", state="
213  << (arg1.state() == IteratorState::valid ?
214  "valid" :
215  (arg1.state() == IteratorState::past_the_end ?
216  "past_the_end" :
217  "invalid")));
222  "Iterators can only be compared if they point to the same "
223  "triangulation, or if neither of them are associated "
224  "with a triangulation.");
225  // TODO: Write documentation!
230  // TODO: Write documentation!
250  // TODO: Write documentation!
256  int,
257  << "You can only set the child index of an even numbered child."
258  << "The number of the child given was " << arg1 << ".");
259 } // namespace TriaAccessorExceptions
260 
261 
288 template <int structdim, int dim, int spacedim = dim>
289 class TriaAccessorBase
290 {
291 public:
297  static const unsigned int space_dimension = spacedim;
298 
304  static const unsigned int dimension = dim;
305 
311  static const unsigned int structure_dimension = structdim;
312 
322  void
323  operator=(const TriaAccessorBase *) = delete;
324 
325 protected:
331  using AccessorData = void;
332 
336  TriaAccessorBase(const Triangulation<dim, spacedim> *parent = nullptr,
337  const int level = -1,
338  const int index = -1,
339  const AccessorData * = nullptr);
340 
345 
353  void
354  copy_from(const TriaAccessorBase &);
355 
360  operator=(const TriaAccessorBase &);
361 
372  bool
373  operator<(const TriaAccessorBase &other) const;
374 
375 protected:
379  bool
380  operator==(const TriaAccessorBase &) const;
381 
385  bool
386  operator!=(const TriaAccessorBase &) const;
387 
401  void
402  operator++();
403 
411  void
412  operator--();
422  objects() const;
423 
424 public:
430  using LocalData = void *;
431 
455  int
456  level() const;
457 
484  int
485  index() const;
486 
492  state() const;
493 
499  get_triangulation() const;
500 
504 protected:
509  typename ::internal::TriaAccessorImplementation::
510  PresentLevelType<structdim, dim>::type present_level;
511 
517 
522 
523 private:
524  template <typename Accessor>
525  friend class TriaRawIterator;
526  template <typename Accessor>
527  friend class TriaIterator;
528  template <typename Accessor>
529  friend class TriaActiveIterator;
530 };
531 
532 
533 
555 template <int structdim, int dim, int spacedim = dim>
556 class InvalidAccessor : public TriaAccessorBase<structdim, dim, spacedim>
557 {
558 public:
562  using AccessorData =
564 
572  InvalidAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
573  const int level = -1,
574  const int index = -1,
575  const AccessorData * local_data = nullptr);
576 
585 
590  template <typename OtherAccessor>
591  InvalidAccessor(const OtherAccessor &);
592 
596  void
597  copy_from(const InvalidAccessor &);
598 
602  bool
603  operator==(const InvalidAccessor &) const;
604  bool
605  operator!=(const InvalidAccessor &) const;
606 
610  void
611  operator++() const;
612  void
613  operator--() const;
614 
619  bool
620  used() const;
621 
626  bool
627  has_children() const;
628 
633  manifold_id() const;
634 
638  unsigned int
639  user_index() const;
640 
644  void
645  set_user_index(const unsigned int p) const;
646 
650  void
651  set_manifold_id(const types::manifold_id) const;
652 
657  vertex(const unsigned int i) const;
658 
663  typename ::internal::TriangulationImplementation::
664  Iterators<dim, spacedim>::line_iterator
665  line(const unsigned int i) const;
666 
671  typename ::internal::TriangulationImplementation::
672  Iterators<dim, spacedim>::quad_iterator
673  quad(const unsigned int i) const;
674 };
675 
676 
677 
696 template <int structdim, int dim, int spacedim>
697 class TriaAccessor : public TriaAccessorBase<structdim, dim, spacedim>
698 {
699 public:
703  using AccessorData =
705 
709  TriaAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
710  const int level = -1,
711  const int index = -1,
712  const AccessorData * local_data = nullptr);
713 
726  template <int structdim2, int dim2, int spacedim2>
728 
733  template <int structdim2, int dim2, int spacedim2>
735 
745  void
746  operator=(const TriaAccessor &) = delete;
747 
754  bool
755  used() const;
756 
769  vertex_iterator(const unsigned int i) const;
770 
786  unsigned int
787  vertex_index(const unsigned int i) const;
788 
827  vertex(const unsigned int i) const;
828 
832  typename ::internal::TriangulationImplementation::
833  Iterators<dim, spacedim>::line_iterator
834  line(const unsigned int i) const;
835 
842  unsigned int
843  line_index(const unsigned int i) const;
844 
848  typename ::internal::TriangulationImplementation::
849  Iterators<dim, spacedim>::quad_iterator
850  quad(const unsigned int i) const;
851 
858  unsigned int
859  quad_index(const unsigned int i) const;
882  bool
883  face_orientation(const unsigned int face) const;
884 
894  bool
895  face_flip(const unsigned int face) const;
896 
906  bool
907  face_rotation(const unsigned int face) const;
908 
919  bool
920  line_orientation(const unsigned int line) const;
935  bool
936  has_children() const;
937 
942  unsigned int
943  n_children() const;
944 
958  unsigned int
959  number_of_children() const;
960 
974  unsigned int
975  max_refinement_depth() const;
976 
981  child(const unsigned int i) const;
982 
992  isotropic_child(const unsigned int i) const;
993 
998  refinement_case() const;
999 
1005  int
1006  child_index(const unsigned int i) const;
1007 
1013  int
1014  isotropic_child_index(const unsigned int i) const;
1037  boundary_id() const;
1038 
1068  void
1069  set_boundary_id(const types::boundary_id) const;
1070 
1101  void
1103 
1111  bool
1112  at_boundary() const;
1113 
1123  const Manifold<dim, spacedim> &
1124  get_manifold() const;
1125 
1147  manifold_id() const;
1148 
1166  void
1167  set_manifold_id(const types::manifold_id) const;
1168 
1182  void
1184 
1201  bool
1202  user_flag_set() const;
1203 
1209  void
1210  set_user_flag() const;
1211 
1217  void
1218  clear_user_flag() const;
1219 
1225  void
1226  recursively_set_user_flag() const;
1227 
1233  void
1235 
1241  void
1242  clear_user_data() const;
1243 
1255  void
1256  set_user_pointer(void *p) const;
1257 
1263  void
1264  clear_user_pointer() const;
1265 
1281  void *
1282  user_pointer() const;
1283 
1305  void
1306  recursively_set_user_pointer(void *p) const;
1307 
1314  void
1316 
1326  void
1327  set_user_index(const unsigned int p) const;
1328 
1334  void
1335  clear_user_index() const;
1336 
1348  unsigned int
1349  user_index() const;
1350 
1368  void
1369  recursively_set_user_index(const unsigned int p) const;
1370 
1379  void
1399  double
1400  diameter() const;
1401 
1428  std::pair<Point<spacedim>, double>
1429  enclosing_ball() const;
1430 
1440  bounding_box() const;
1441 
1451  double
1452  extent_in_direction(const unsigned int axis) const;
1453 
1457  double
1458  minimum_vertex_distance() const;
1459 
1474  intermediate_point(const Point<structdim> &coordinates) const;
1475 
1500 
1536  center(const bool respect_manifold = false,
1537  const bool interpolate_from_surrounding = false) const;
1538 
1557  barycenter() const;
1558 
1584  double
1585  measure() const;
1586 
1601  bool
1604 
1610 private:
1615  void
1617 
1622  void
1623  set(const ::internal::TriangulationImplementation::TriaObject<structdim>
1624  &o) const;
1625 
1633  void
1634  set_line_orientation(const unsigned int line, const bool orientation) const;
1635 
1646  void
1647  set_face_orientation(const unsigned int face, const bool orientation) const;
1648 
1655  void
1656  set_face_flip(const unsigned int face, const bool flip) const;
1657 
1664  void
1665  set_face_rotation(const unsigned int face, const bool rotation) const;
1666 
1670  void
1671  set_used_flag() const;
1672 
1676  void
1677  clear_used_flag() const;
1678 
1687  void
1688  set_refinement_case(const RefinementCase<structdim> &ref_case) const;
1689 
1697  void
1698  clear_refinement_case() const;
1699 
1706  void
1707  set_children(const unsigned int i, const int index) const;
1708 
1713  void
1714  clear_children() const;
1715 
1716 private:
1717  template <int, int>
1718  friend class Triangulation;
1719 
1720  friend struct ::internal::TriangulationImplementation::Implementation;
1721  friend struct ::internal::TriaAccessorImplementation::Implementation;
1722 };
1723 
1724 
1725 
1745 template <int dim, int spacedim>
1746 class TriaAccessor<0, dim, spacedim>
1747 {
1748 public:
1754  static const unsigned int space_dimension = spacedim;
1755 
1761  static const unsigned int dimension = dim;
1762 
1768  static const unsigned int structure_dimension = 0;
1769 
1773  using AccessorData = void;
1774 
1780  const unsigned int vertex_index);
1781 
1788  const int level = 0,
1789  const int index = 0,
1790  const AccessorData * = nullptr);
1791 
1795  template <int structdim2, int dim2, int spacedim2>
1797 
1801  template <int structdim2, int dim2, int spacedim2>
1803 
1808  state() const;
1809 
1814  static int
1815  level();
1816 
1821  int
1822  index() const;
1823 
1829  get_triangulation() const;
1830 
1840  void
1841  operator++();
1842 
1846  void
1847  operator--();
1851  bool
1852  operator==(const TriaAccessor &) const;
1853 
1857  bool
1858  operator!=(const TriaAccessor &) const;
1859 
1887  unsigned int
1888  vertex_index(const unsigned int i = 0) const;
1889 
1895  Point<spacedim> &
1896  vertex(const unsigned int i = 0) const;
1897 
1902  typename ::internal::TriangulationImplementation::
1903  Iterators<dim, spacedim>::line_iterator static line(const unsigned int);
1904 
1908  static unsigned int
1909  line_index(const unsigned int i);
1910 
1914  static typename ::internal::TriangulationImplementation::
1915  Iterators<dim, spacedim>::quad_iterator
1916  quad(const unsigned int i);
1917 
1921  static unsigned int
1922  quad_index(const unsigned int i);
1923 
1939  double
1940  diameter() const;
1941 
1949  double
1950  extent_in_direction(const unsigned int axis) const;
1951 
1960  center(const bool respect_manifold = false,
1961  const bool interpolate_from_surrounding = false) const;
1962 
1971  double
1972  measure() const;
1987  static bool
1988  face_orientation(const unsigned int face);
1989 
1993  static bool
1994  face_flip(const unsigned int face);
1995 
1999  static bool
2000  face_rotation(const unsigned int face);
2001 
2005  static bool
2006  line_orientation(const unsigned int line);
2007 
2022  static bool
2023  has_children();
2024 
2029  static unsigned int
2030  n_children();
2031 
2036  static unsigned int
2038 
2042  static unsigned int
2044 
2049  child(const unsigned int);
2050 
2055  isotropic_child(const unsigned int);
2056 
2060  static RefinementCase<0>
2061  refinement_case();
2062 
2066  static int
2067  child_index(const unsigned int i);
2068 
2072  static int
2073  isotropic_child_index(const unsigned int i);
2081  bool
2082  used() const;
2083 
2084 protected:
2092  void
2093  copy_from(const TriaAccessor &);
2094 
2103  bool
2104  operator<(const TriaAccessor &other) const;
2105 
2110 
2114  unsigned int global_vertex_index;
2115 
2116 private:
2117  template <typename Accessor>
2118  friend class TriaRawIterator;
2119  template <typename Accessor>
2120  friend class TriaIterator;
2121  template <typename Accessor>
2122  friend class TriaActiveIterator;
2123 };
2124 
2125 
2126 
2144 template <int spacedim>
2145 class TriaAccessor<0, 1, spacedim>
2146 {
2147 public:
2153  static const unsigned int space_dimension = spacedim;
2154 
2160  static const unsigned int dimension = 1;
2161 
2167  static const unsigned int structure_dimension = 0;
2168 
2172  using AccessorData = void;
2173 
2179  {
2191  right_vertex
2192  };
2193 
2206  const VertexKind vertex_kind,
2207  const unsigned int vertex_index);
2208 
2214  TriaAccessor(const Triangulation<1, spacedim> *tria = nullptr,
2215  const int = 0,
2216  const int = 0,
2217  const AccessorData * = nullptr);
2218 
2222  template <int structdim2, int dim2, int spacedim2>
2224 
2228  template <int structdim2, int dim2, int spacedim2>
2230 
2235  void
2236  copy_from(const TriaAccessor &);
2237 
2244  state();
2245 
2250  static int
2251  level();
2252 
2257  int
2258  index() const;
2259 
2265  get_triangulation() const;
2266 
2277  void
2278  operator++() const;
2279 
2284  void
2285  operator--() const;
2289  bool
2290  operator==(const TriaAccessor &) const;
2291 
2295  bool
2296  operator!=(const TriaAccessor &) const;
2297 
2306  bool
2307  operator<(const TriaAccessor &other) const;
2308 
2335  unsigned int
2336  vertex_index(const unsigned int i = 0) const;
2337 
2343  Point<spacedim> &
2344  vertex(const unsigned int i = 0) const;
2345 
2351  center() const;
2352 
2357  typename ::internal::TriangulationImplementation::
2358  Iterators<1, spacedim>::line_iterator static line(const unsigned int);
2359 
2366  static unsigned int
2367  line_index(const unsigned int i);
2368 
2372  static typename ::internal::TriangulationImplementation::
2373  Iterators<1, spacedim>::quad_iterator
2374  quad(const unsigned int i);
2375 
2382  static unsigned int
2383  quad_index(const unsigned int i);
2384 
2394  bool
2395  at_boundary() const;
2396 
2412  boundary_id() const;
2413 
2417  const Manifold<1, spacedim> &
2418  get_manifold() const;
2419 
2427  manifold_id() const;
2428 
2429 
2440  static bool
2441  face_orientation(const unsigned int face);
2442 
2446  static bool
2447  face_flip(const unsigned int face);
2448 
2452  static bool
2453  face_rotation(const unsigned int face);
2454 
2458  static bool
2459  line_orientation(const unsigned int line);
2460 
2475  static bool
2476  has_children();
2477 
2482  static unsigned int
2483  n_children();
2484 
2489  static unsigned int
2491 
2495  static unsigned int
2497 
2502  child(const unsigned int);
2503 
2508  isotropic_child(const unsigned int);
2509 
2513  static RefinementCase<0>
2514  refinement_case();
2515 
2519  static int
2520  child_index(const unsigned int i);
2521 
2525  static int
2526  isotropic_child_index(const unsigned int i);
2559  void
2561 
2568  void
2570 
2582  void
2584 
2596  void
2605  bool
2606  used() const;
2607 
2608 protected:
2613 
2619 
2623  unsigned int global_vertex_index;
2624 };
2625 
2626 
2627 
2644 template <int dim, int spacedim = dim>
2645 class CellAccessor : public TriaAccessor<dim, dim, spacedim>
2646 {
2647 public:
2652 
2657 
2669  const int level = -1,
2670  const int index = -1,
2671  const AccessorData * local_data = nullptr);
2672 
2676  CellAccessor(const TriaAccessor<dim, dim, spacedim> &cell_accessor);
2677 
2690  template <int structdim2, int dim2, int spacedim2>
2692 
2697  template <int structdim2, int dim2, int spacedim2>
2699 
2709  void
2710  operator=(const CellAccessor<dim, spacedim> &) = delete;
2711 
2728  child(const unsigned int i) const;
2729 
2733  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>
2734  face(const unsigned int i) const;
2735 
2745  unsigned int
2746  face_index(const unsigned int i) const;
2747 
2796  neighbor_child_on_subface(const unsigned int face_no,
2797  const unsigned int subface_no) const;
2798 
2821  neighbor(const unsigned int i) const;
2822 
2827  int
2828  neighbor_index(const unsigned int i) const;
2829 
2834  int
2835  neighbor_level(const unsigned int i) const;
2836 
2848  unsigned int
2849  neighbor_of_neighbor(const unsigned int neighbor) const;
2850 
2861  bool
2862  neighbor_is_coarser(const unsigned int neighbor) const;
2863 
2878  std::pair<unsigned int, unsigned int>
2879  neighbor_of_coarser_neighbor(const unsigned int neighbor) const;
2880 
2887  unsigned int
2888  neighbor_face_no(const unsigned int neighbor) const;
2889 
2893  static bool
2894  is_level_cell();
2895 
2909  bool
2910  has_periodic_neighbor(const unsigned int i) const;
2911 
2929  periodic_neighbor(const unsigned int i) const;
2930 
2939  neighbor_or_periodic_neighbor(const unsigned int i) const;
2940 
2956  periodic_neighbor_child_on_subface(const unsigned int face_no,
2957  const unsigned int subface_no) const;
2958 
2969  std::pair<unsigned int, unsigned int>
2970  periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const;
2971 
2977  int
2978  periodic_neighbor_index(const unsigned int i) const;
2979 
2985  int
2986  periodic_neighbor_level(const unsigned int i) const;
2987 
3002  unsigned int
3003  periodic_neighbor_of_periodic_neighbor(const unsigned int i) const;
3004 
3010  unsigned int
3011  periodic_neighbor_face_no(const unsigned int i) const;
3012 
3019  bool
3020  periodic_neighbor_is_coarser(const unsigned int i) const;
3021 
3038  bool
3039  at_boundary(const unsigned int i) const;
3040 
3049  bool
3050  at_boundary() const;
3051 
3059  bool
3060  has_boundary_lines() const;
3087  refine_flag_set() const;
3088 
3106  void
3107  set_refine_flag(const RefinementCase<dim> ref_case =
3109 
3113  void
3114  clear_refine_flag() const;
3115 
3123  bool
3125  const unsigned int face_no,
3126  const RefinementCase<dim - 1> &face_refinement_case =
3128 
3134  bool
3135  flag_for_line_refinement(const unsigned int line_no) const;
3136 
3146  subface_case(const unsigned int face_no) const;
3147 
3151  bool
3152  coarsen_flag_set() const;
3153 
3158  void
3159  set_coarsen_flag() const;
3160 
3164  void
3165  clear_coarsen_flag() const;
3189  material_id() const;
3190 
3202  void
3203  set_material_id(const types::material_id new_material_id) const;
3204 
3213  void
3214  recursively_set_material_id(const types::material_id new_material_id) const;
3241  subdomain_id() const;
3242 
3258  void
3259  set_subdomain_id(const types::subdomain_id new_subdomain_id) const;
3260 
3266  level_subdomain_id() const;
3267 
3272  void
3274  const types::subdomain_id new_level_subdomain_id) const;
3275 
3276 
3292  void
3294  const types::subdomain_id new_subdomain_id) const;
3312  bool
3313  direction_flag() const;
3314 
3340  unsigned int
3341  active_cell_index() const;
3342 
3350  int
3351  parent_index() const;
3352 
3359  parent() const;
3360 
3380  bool
3381  active() const;
3382 
3402  bool
3403  is_locally_owned() const;
3404 
3409  bool
3410  is_locally_owned_on_level() const;
3411 
3435  bool
3436  is_ghost() const;
3437 
3464  bool
3465  is_artificial() const;
3466 
3480  bool
3481  point_inside(const Point<spacedim> &p) const;
3482 
3491  void
3492  set_neighbor(const unsigned int i,
3493  const TriaIterator<CellAccessor<dim, spacedim>> &pointer) const;
3494 
3508  CellId
3509  id() const;
3510 
3528 
3529 protected:
3545  unsigned int
3546  neighbor_of_neighbor_internal(const unsigned int neighbor) const;
3547 
3553  template <int dim_, int spacedim_>
3554  bool
3555  point_inside_codim(const Point<spacedim_> &p) const;
3556 
3557 
3558 
3559 private:
3564  void
3565  set_active_cell_index(const unsigned int active_cell_index);
3566 
3570  void
3571  set_parent(const unsigned int parent_index);
3572 
3579  void
3580  set_direction_flag(const bool new_direction_flag) const;
3581 
3582  template <int, int>
3583  friend class Triangulation;
3584 
3585  friend struct ::internal::TriangulationImplementation::Implementation;
3586 };
3587 
3588 
3589 
3590 /* ----- declaration of explicit specializations and general templates ----- */
3591 
3592 
3593 template <int structdim, int dim, int spacedim>
3594 template <typename OtherAccessor>
3596  const OtherAccessor &)
3597 {
3598  Assert(false,
3599  ExcMessage("You are attempting an illegal conversion between "
3600  "iterator/accessor types. The constructor you call "
3601  "only exists to make certain template constructs "
3602  "easier to write as dimension independent code but "
3603  "the conversion is not valid in the current context."));
3604 }
3605 
3606 
3607 
3608 template <int structdim, int dim, int spacedim>
3609 template <int structdim2, int dim2, int spacedim2>
3612 {
3613  Assert(false,
3614  ExcMessage("You are attempting an illegal conversion between "
3615  "iterator/accessor types. The constructor you call "
3616  "only exists to make certain template constructs "
3617  "easier to write as dimension independent code but "
3618  "the conversion is not valid in the current context."));
3619 }
3620 
3621 
3622 
3623 template <int dim, int spacedim>
3624 template <int structdim2, int dim2, int spacedim2>
3627 {
3628  Assert(false,
3629  ExcMessage("You are attempting an illegal conversion between "
3630  "iterator/accessor types. The constructor you call "
3631  "only exists to make certain template constructs "
3632  "easier to write as dimension independent code but "
3633  "the conversion is not valid in the current context."));
3634 }
3635 
3636 
3637 
3638 template <int structdim, int dim, int spacedim>
3639 template <int structdim2, int dim2, int spacedim2>
3642 {
3643  Assert(false,
3644  ExcMessage("You are attempting an illegal conversion between "
3645  "iterator/accessor types. The constructor you call "
3646  "only exists to make certain template constructs "
3647  "easier to write as dimension independent code but "
3648  "the conversion is not valid in the current context."));
3649 }
3650 
3651 
3652 
3653 template <int dim, int spacedim>
3654 template <int structdim2, int dim2, int spacedim2>
3657 {
3658  Assert(false,
3659  ExcMessage("You are attempting an illegal conversion between "
3660  "iterator/accessor types. The constructor you call "
3661  "only exists to make certain template constructs "
3662  "easier to write as dimension independent code but "
3663  "the conversion is not valid in the current context."));
3664 }
3665 
3666 
3667 #ifndef DOXYGEN
3668 
3669 template <>
3670 bool
3672 template <>
3673 bool
3675 template <>
3676 bool
3678 template <>
3679 bool
3681 template <>
3682 bool
3684 template <>
3685 bool
3687 // -------------------------------------------------------------------
3688 
3689 template <>
3690 void
3692 
3693 #endif // DOXYGEN
3694 
3695 DEAL_II_NAMESPACE_CLOSE
3696 
3697 // include more templates in debug and optimized mode
3698 #include "tria_accessor.templates.h"
3699 
3700 #endif
unsigned int periodic_neighbor_of_periodic_neighbor(const unsigned int i) const
void clear_used_flag() const
std::pair< Point< spacedim >, double > enclosing_ball() const
unsigned int manifold_id
Definition: types.h:123
static ::ExceptionBase & ExcCellHasNoChildren()
InvalidAccessor(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
bool point_inside_codim(const Point< spacedim_ > &p) const
void set_face_rotation(const unsigned int face, const bool rotation) const
TriaIterator< CellAccessor< dim, spacedim > > neighbor_or_periodic_neighbor(const unsigned int i) const
static ::ExceptionBase & ExcNeighborIsNotCoarser()
const Triangulation< dim, spacedim > * tria
void set_active_cell_index(const unsigned int active_cell_index)
void recursively_clear_user_pointer() const
unsigned int active_cell_index() const
typename ::internal::TriangulationImplementation::Iterators< dim, spacedim >::line_iterator line(const unsigned int i) const
void set_manifold_id(const types::manifold_id) const
void set_neighbor(const unsigned int i, const TriaIterator< CellAccessor< dim, spacedim >> &pointer) const
bool operator==(const InvalidAccessor &) const
static ::ExceptionBase & ExcNeighborIsCoarser()
static ::ExceptionBase & ExcRefineCellNotActive()
unsigned int neighbor_of_neighbor_internal(const unsigned int neighbor) const
unsigned int user_index() const
static ::ExceptionBase & ExcDereferenceInvalidObject(AccessorType arg1)
unsigned int user_index() const
static bool is_level_cell()
int isotropic_child_index(const unsigned int i) const
void set_user_index(const unsigned int p) const
const Manifold< dim, spacedim > & get_manifold() const
void * user_pointer() const
types::material_id material_id() const
bool at_boundary() const
bool flag_for_line_refinement(const unsigned int line_no) const
void set_line_orientation(const unsigned int line, const bool orientation) const
TriaIterator< TriaAccessor< structdim, dim, spacedim > > isotropic_child(const unsigned int i) const
void copy_from(const TriaAccessorBase &)
double extent_in_direction(const unsigned int axis) const
const Triangulation< 1, spacedim > * tria
bool operator!=(const TriaAccessorBase &) const
unsigned int material_id
Definition: types.h:134
unsigned int n_children() const
void set_all_manifold_ids(const types::manifold_id) const
bool is_translation_of(const TriaIterator< TriaAccessor< structdim, dim, spacedim >> &o) const
void set_children(const unsigned int i, const int index) const
void clear_user_pointer() const
unsigned int line_index(const unsigned int i) const
static ::ExceptionBase & ExcFacesHaveNoLevel()
void operator=(const CellAccessor< dim, spacedim > &)=delete
unsigned int face_index(const unsigned int i) const
TriaIterator< CellAccessor< dim, spacedim > > neighbor_child_on_subface(const unsigned int face_no, const unsigned int subface_no) const
typename TriaAccessorBase< structdim, dim, spacedim >::AccessorData AccessorData
bool neighbor_is_coarser(const unsigned int neighbor) const
bool has_children() const
void operator++() const
int level() const
unsigned int vertex_index(const unsigned int i) const
types::manifold_id manifold_id() const
bool is_locally_owned() const
bool face_rotation(const unsigned int face) const
static ::ExceptionBase & ExcCellNotUsed()
void copy_from(const InvalidAccessor &)
bool at_boundary() const
void clear_user_data() const
TriaIterator< CellAccessor< dim, spacedim > > parent() const
static const unsigned int space_dimension
CellId id() const
void set_user_pointer(void *p) const
void set_boundary_id_internal(const types::boundary_id id) const
bool coarsen_flag_set() const
TriaIterator< CellAccessor< dim, spacedim > > neighbor(const unsigned int i) const
typename ::internal::TriangulationImplementation::Iterators< dim, spacedim >::quad_iterator quad(const unsigned int i) const
static ::ExceptionBase & ExcMessage(std::string arg1)
bool has_periodic_neighbor(const unsigned int i) const
void recursively_set_user_flag() const
bool operator==(const TriaAccessorBase &) const
void set_coarsen_flag() const
unsigned int subdomain_id
Definition: types.h:43
TriaIterator< TriaAccessor< dim - 1, dim, spacedim > > face(const unsigned int i) const
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:518
void clear_user_flag() const
types::manifold_id manifold_id() const
TriaIterator< TriaAccessor< 0, dim, spacedim > > vertex_iterator(const unsigned int i) const
const Triangulation< dim, spacedim > * tria
bool active() const
bool is_locally_owned_on_level() const
bool periodic_neighbor_is_coarser(const unsigned int i) const
static ::ExceptionBase & ExcCellFlaggedForCoarsening()
static ::ExceptionBase & ExcSetOnlyEvenChildren(int arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1407
std::pair< unsigned int, unsigned int > periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const
Point< spacedim > & vertex(const unsigned int i) const
BoundingBox< spacedim > bounding_box() const
bool line_orientation(const unsigned int line) const
void set_material_id(const types::material_id new_material_id) const
bool is_ghost() const
Point< spacedim > & vertex(const unsigned int i) const
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:496
#define DeclException0(Exception0)
Definition: exceptions.h:473
void set_face_orientation(const unsigned int face, const bool orientation) const
unsigned int quad_index(const unsigned int i) const
bool operator<(const TriaAccessorBase &other) const
RefinementCase< dim > refine_flag_set() const
int index() const
static ::ExceptionBase & ExcCellNotActive()
void operator=(const TriaAccessorBase *)=delete
void operator=(const TriaAccessor &)=delete
IteratorState::IteratorStates state() const
void set_refine_flag(const RefinementCase< dim > ref_case=RefinementCase< dim >::isotropic_refinement) const
static ::ExceptionBase & ExcNoPeriodicNeighbor()
static ::ExceptionBase & ExcCantSetChildren(int arg1)
void clear_refine_flag() const
void recursively_set_material_id(const types::material_id new_material_id) const
void clear_children() const
void recursively_clear_user_index() const
TriaIterator< CellAccessor< dim, spacedim > > periodic_neighbor_child_on_subface(const unsigned int face_no, const unsigned int subface_no) const
void set_user_index(const unsigned int p) const
double minimum_vertex_distance() const
void clear_user_index() const
bool has_children() const
void recursively_clear_user_flag() const
Definition: cell_id.h:63
CellAccessor(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
Point< spacedim > center(const bool respect_manifold=false, const bool interpolate_from_surrounding=false) const
int parent_index() const
TriaAccessorBase(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *=nullptr)
static ::ExceptionBase & ExcCellHasNoParent()
void set_all_boundary_ids(const types::boundary_id) const
typename ::internal::TriangulationImplementation::Iterators< dim, spacedim >::quad_iterator quad(const unsigned int i) const
Point< structdim > real_to_unit_cell_affine_approximation(const Point< spacedim > &point) const
bool is_artificial() const
types::boundary_id boundary_id() const
double diameter() const
int periodic_neighbor_level(const unsigned int i) const
const Triangulation< dim, spacedim > & get_triangulation() const
void recursively_set_user_pointer(void *p) const
bool flag_for_face_refinement(const unsigned int face_no, const RefinementCase< dim - 1 > &face_refinement_case=RefinementCase< dim - 1 >::isotropic_refinement) const
std::pair< unsigned int, unsigned int > neighbor_of_coarser_neighbor(const unsigned int neighbor) const
void set_manifold_id(const types::manifold_id) const
void recursively_set_subdomain_id(const types::subdomain_id new_subdomain_id) const
bool face_flip(const unsigned int face) const
::internal::SubfaceCase< dim > subface_case(const unsigned int face_no) const
int neighbor_level(const unsigned int i) const
bool point_inside(const Point< spacedim > &p) const
Point< spacedim > barycenter() const
TriaIterator< TriaAccessor< structdim, dim, spacedim > > child(const unsigned int i) const
bool used() const
void set_direction_flag(const bool new_direction_flag) const
static const unsigned int structure_dimension
unsigned int periodic_neighbor_face_no(const unsigned int i) const
int periodic_neighbor_index(const unsigned int i) const
TriaIterator< CellAccessor< dim, spacedim > > periodic_neighbor(const unsigned int i) const
unsigned int number_of_children() const
void clear_refinement_case() const
int child_index(const unsigned int i) const
void set_parent(const unsigned int parent_index)
void clear_coarsen_flag() const
bool used() const
bool face_orientation(const unsigned int face) const
Point< spacedim > intermediate_point(const Point< structdim > &coordinates) const
void set_boundary_id(const types::boundary_id) const
TriaAccessor(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
bool direction_flag() const
int neighbor_index(const unsigned int i) const
Iterator reached end of container.
Iterator points to a valid object.
void set_level_subdomain_id(const types::subdomain_id new_level_subdomain_id) const
void recursively_set_user_index(const unsigned int p) const
void set_subdomain_id(const types::subdomain_id new_subdomain_id) const
types::subdomain_id level_subdomain_id() const
void set_user_flag() const
bool user_flag_set() const
void set_face_flip(const unsigned int face, const bool flip) const
unsigned int neighbor_face_no(const unsigned int neighbor) const
types::subdomain_id subdomain_id() const
static ::ExceptionBase & ExcCellFlaggedForRefinement()
static const unsigned int dimension
double measure() const
void set_refinement_case(const RefinementCase< structdim > &ref_case) const
RefinementCase< structdim > refinement_case() const
unsigned int max_refinement_depth() const
unsigned int boundary_id
Definition: types.h:111
unsigned int neighbor_of_neighbor(const unsigned int neighbor) const
typename ::internal::TriaAccessorImplementation::PresentLevelType< structdim, dim >::type present_level
::internal::TriangulationImplementation::TriaObjects< ::internal::TriangulationImplementation::TriaObject< structdim > > & objects() const
void set_used_flag() const
TriaIterator< CellAccessor< dim, spacedim > > child(const unsigned int i) const
bool has_boundary_lines() const
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcCantCompareIterators()
typename ::internal::TriangulationImplementation::Iterators< dim, spacedim >::line_iterator line(const unsigned int i) const