Reference documentation for deal.II version Git e159ac89de 2020-06-06 19:38:41 +0200
\(\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 - 2020 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>
30 
31 #include <utility>
32 
33 
35 
36 // Forward declarations
37 #ifndef DOXYGEN
38 template <int dim, int spacedim>
39 class Triangulation;
40 template <typename Accessor>
41 class TriaRawIterator;
42 template <typename Accessor>
43 class TriaIterator;
44 template <typename Accessor>
45 class TriaActiveIterator;
46 
47 template <int dim, int spacedim>
48 class Manifold;
49 #endif
50 
51 namespace internal
52 {
53  namespace TriangulationImplementation
54  {
55  class TriaObjects;
56  struct Implementation;
57  } // namespace TriangulationImplementation
58 
59  namespace TriaAccessorImplementation
60  {
61  struct Implementation;
62 
68  template <int structdim, int dim>
70  {
71  struct type
72  {
76  type() = default;
77 
81  type(const int level)
82  {
83  Assert(level == 0, ExcInternalError());
84  (void)level; // removes -Wunused-parameter warning in optimized mode
85  }
86 
90  operator int() const
91  {
92  return 0;
93  }
94 
95  void
96  operator++() const
97  {
98  Assert(false, ExcInternalError());
99  }
100 
101  void
102  operator--() const
103  {
104  Assert(false, ExcInternalError());
105  }
106  };
107  };
108 
109 
115  template <int dim>
116  struct PresentLevelType<dim, dim>
117  {
118  using type = int;
119  };
120 
121  } // namespace TriaAccessorImplementation
122 } // namespace internal
123 template <int structdim, int dim, int spacedim>
125 template <int dim, int spacedim>
126 class TriaAccessor<0, dim, spacedim>;
127 template <int spacedim>
128 class TriaAccessor<0, 1, spacedim>;
129 
130 // note: the file tria_accessor.templates.h is included at the end of
131 // this file. this includes a lot of templates. originally, this was
132 // only done in debug mode, but led to cyclic reduction problems and
133 // so is now on by default.
134 
135 
140 {
145  "The operation you are attempting can only be performed for "
146  "(cell, face, or edge) iterators that point to valid "
147  "objects. These objects need not necessarily be active, "
148  "i.e., have no children, but they need to be part of a "
149  "triangulation. (The objects pointed to by an iterator "
150  "may -- after coarsening -- also be objects that used "
151  "to be part of a triangulation, but are now no longer "
152  "used. Their memory location may have been retained "
153  "for re-use upon the next mesh refinement, but is "
154  "currently unused.)");
165  "The operation you are attempting can only be performed for "
166  "(cell, face, or edge) iterators that point to 'active' "
167  "objects. 'Active' objects are those that do not have "
168  "children (in the case of cells), or that are part of "
169  "an active cell (in the case of faces or edges). However, "
170  "the object on which you are trying the current "
171  "operation is not 'active' in this sense.");
178  "The operation you are attempting can only be performed for "
179  "(cell, face, or edge) iterators that have children, "
180  "but the object on which you are trying the current "
181  "operation does not have any.");
189  "The operation you are attempting can only be performed for "
190  "(cell, face, or edge) iterators that have a parent object, "
191  "but the object on which you are trying the current "
192  "operation does not have one -- i.e., it is on the "
193  "coarsest level of the triangulation.");
198  int,
199  << "You can only set the child index if the cell does not "
200  << "currently have children registered; or you can clear it. "
201  << "The given index was " << arg1
202  << " (-1 means: clear children).");
206  template <typename AccessorType>
208  AccessorType,
209  << "You tried to dereference an iterator for which this "
210  << "is not possible. More information on this iterator: "
211  << "index=" << arg1.index() << ", state="
212  << (arg1.state() == IteratorState::valid ?
213  "valid" :
214  (arg1.state() == IteratorState::past_the_end ?
215  "past_the_end" :
216  "invalid")));
221  "Iterators can only be compared if they point to the same "
222  "triangulation, or if neither of them are associated "
223  "with a triangulation.");
224  // TODO: Write documentation!
229  // TODO: Write documentation!
249  // TODO: Write documentation!
255  int,
256  << "You can only set the child index of an even numbered child."
257  << "The number of the child given was " << arg1 << ".");
258 } // namespace TriaAccessorExceptions
259 
260 
286 template <int structdim, int dim, int spacedim = dim>
288 {
289 public:
295  static const unsigned int space_dimension = spacedim;
296 
302  static const unsigned int dimension = dim;
303 
309  static const unsigned int structure_dimension = structdim;
310 
320  void
321  operator=(const TriaAccessorBase *) = delete;
322 
323 protected:
329  using AccessorData = void;
330 
334  TriaAccessorBase(const Triangulation<dim, spacedim> *parent = nullptr,
335  const int level = -1,
336  const int index = -1,
337  const AccessorData * = nullptr);
338 
343 
351  void
352  copy_from(const TriaAccessorBase &);
353 
358  operator=(const TriaAccessorBase &);
359 
370  bool
371  operator<(const TriaAccessorBase &other) const;
372 
373 protected:
377  bool
378  operator==(const TriaAccessorBase &) const;
379 
383  bool
384  operator!=(const TriaAccessorBase &) const;
385 
399  void
400  operator++();
401 
409  void
410  operator--();
419  objects() const;
420 
421 public:
427  using LocalData = void *;
428 
452  int
453  level() const;
454 
481  int
482  index() const;
483 
489  state() const;
490 
496  get_triangulation() const;
497 
501 protected:
506  typename ::internal::TriaAccessorImplementation::
507  PresentLevelType<structdim, dim>::type present_level;
508 
514 
519 
520 private:
521  template <typename Accessor>
522  friend class TriaRawIterator;
523  template <typename Accessor>
524  friend class TriaIterator;
525  template <typename Accessor>
526  friend class TriaActiveIterator;
527 };
528 
529 
530 
551 template <int structdim, int dim, int spacedim = dim>
552 class InvalidAccessor : public TriaAccessorBase<structdim, dim, spacedim>
553 {
554 public:
558  using AccessorData =
560 
568  InvalidAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
569  const int level = -1,
570  const int index = -1,
571  const AccessorData * local_data = nullptr);
572 
581 
586  template <typename OtherAccessor>
587  InvalidAccessor(const OtherAccessor &);
588 
592  void
593  copy_from(const InvalidAccessor &);
594 
598  bool
599  operator==(const InvalidAccessor &) const;
600  bool
601  operator!=(const InvalidAccessor &) const;
602 
606  void
607  operator++() const;
608  void
609  operator--() const;
610 
615  bool
616  used() const;
617 
622  bool
623  has_children() const;
624 
629  manifold_id() const;
630 
634  unsigned int
635  user_index() const;
636 
640  void
641  set_user_index(const unsigned int p) const;
642 
646  void
647  set_manifold_id(const types::manifold_id) const;
648 
653  vertex(const unsigned int i) const;
654 
659  typename ::internal::TriangulationImplementation::
660  Iterators<dim, spacedim>::line_iterator
661  line(const unsigned int i) const;
662 
667  typename ::internal::TriangulationImplementation::
668  Iterators<dim, spacedim>::quad_iterator
669  quad(const unsigned int i) const;
670 };
671 
672 
673 
691 template <int structdim, int dim, int spacedim>
692 class TriaAccessor : public TriaAccessorBase<structdim, dim, spacedim>
693 {
694 public:
698  using AccessorData =
700 
704  TriaAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
705  const int level = -1,
706  const int index = -1,
707  const AccessorData * local_data = nullptr);
708 
721  template <int structdim2, int dim2, int spacedim2>
723 
728  template <int structdim2, int dim2, int spacedim2>
730 
740  void
741  operator=(const TriaAccessor &) = delete;
742 
749  bool
750  used() const;
751 
764  vertex_iterator(const unsigned int i) const;
765 
781  unsigned int
782  vertex_index(const unsigned int i) const;
783 
822  vertex(const unsigned int i) const;
823 
827  typename ::internal::TriangulationImplementation::
828  Iterators<dim, spacedim>::line_iterator
829  line(const unsigned int i) const;
830 
837  unsigned int
838  line_index(const unsigned int i) const;
839 
843  typename ::internal::TriangulationImplementation::
844  Iterators<dim, spacedim>::quad_iterator
845  quad(const unsigned int i) const;
846 
853  unsigned int
854  quad_index(const unsigned int i) const;
877  bool
878  face_orientation(const unsigned int face) const;
879 
889  bool
890  face_flip(const unsigned int face) const;
891 
901  bool
902  face_rotation(const unsigned int face) const;
903 
914  bool
915  line_orientation(const unsigned int line) const;
930  bool
931  has_children() const;
932 
937  unsigned int
938  n_children() const;
939 
953  unsigned int
954  number_of_children() const;
955 
969  unsigned int
970  max_refinement_depth() const;
971 
976  child(const unsigned int i) const;
977 
982  unsigned int
983  child_iterator_to_index(
985 
995  isotropic_child(const unsigned int i) const;
996 
1001  refinement_case() const;
1002 
1008  int
1009  child_index(const unsigned int i) const;
1010 
1016  int
1017  isotropic_child_index(const unsigned int i) const;
1040  boundary_id() const;
1041 
1071  void
1072  set_boundary_id(const types::boundary_id) const;
1073 
1104  void
1105  set_all_boundary_ids(const types::boundary_id) const;
1106 
1114  bool
1115  at_boundary() const;
1116 
1126  const Manifold<dim, spacedim> &
1127  get_manifold() const;
1128 
1150  manifold_id() const;
1151 
1169  void
1170  set_manifold_id(const types::manifold_id) const;
1171 
1185  void
1186  set_all_manifold_ids(const types::manifold_id) const;
1187 
1204  bool
1205  user_flag_set() const;
1206 
1212  void
1213  set_user_flag() const;
1214 
1220  void
1221  clear_user_flag() const;
1222 
1228  void
1229  recursively_set_user_flag() const;
1230 
1236  void
1237  recursively_clear_user_flag() const;
1238 
1244  void
1245  clear_user_data() const;
1246 
1258  void
1259  set_user_pointer(void *p) const;
1260 
1266  void
1267  clear_user_pointer() const;
1268 
1284  void *
1285  user_pointer() const;
1286 
1308  void
1309  recursively_set_user_pointer(void *p) const;
1310 
1317  void
1318  recursively_clear_user_pointer() const;
1319 
1329  void
1330  set_user_index(const unsigned int p) const;
1331 
1337  void
1338  clear_user_index() const;
1339 
1351  unsigned int
1352  user_index() const;
1353 
1371  void
1372  recursively_set_user_index(const unsigned int p) const;
1373 
1382  void
1383  recursively_clear_user_index() const;
1402  double
1403  diameter() const;
1404 
1431  std::pair<Point<spacedim>, double>
1432  enclosing_ball() const;
1433 
1443  bounding_box() const;
1444 
1454  double
1455  extent_in_direction(const unsigned int axis) const;
1456 
1460  double
1461  minimum_vertex_distance() const;
1462 
1477  intermediate_point(const Point<structdim> &coordinates) const;
1478 
1502  real_to_unit_cell_affine_approximation(const Point<spacedim> &point) const;
1503 
1539  center(const bool respect_manifold = false,
1540  const bool interpolate_from_surrounding = false) const;
1541 
1560  barycenter() const;
1561 
1587  double
1588  measure() const;
1589 
1604  bool
1605  is_translation_of(
1607 
1613 private:
1618  void
1619  set_boundary_id_internal(const types::boundary_id id) const;
1620 
1628  void
1629  set_bounding_object_indices(
1630  const std::initializer_list<int> &new_indices) const;
1631 
1635  void
1636  set_bounding_object_indices(
1637  const std::initializer_list<unsigned int> &new_indices) const;
1638 
1646  void
1647  set_line_orientation(const unsigned int line, const bool orientation) const;
1648 
1659  void
1660  set_face_orientation(const unsigned int face, const bool orientation) const;
1661 
1668  void
1669  set_face_flip(const unsigned int face, const bool flip) const;
1670 
1677  void
1678  set_face_rotation(const unsigned int face, const bool rotation) const;
1679 
1683  void
1684  set_used_flag() const;
1685 
1689  void
1690  clear_used_flag() const;
1691 
1700  void
1701  set_refinement_case(const RefinementCase<structdim> &ref_case) const;
1702 
1710  void
1711  clear_refinement_case() const;
1712 
1719  void
1720  set_children(const unsigned int i, const int index) const;
1721 
1726  void
1727  clear_children() const;
1728 
1729 private:
1730  template <int, int>
1731  friend class Triangulation;
1732 
1733  friend struct ::internal::TriangulationImplementation::Implementation;
1734  friend struct ::internal::TriaAccessorImplementation::Implementation;
1735 };
1736 
1737 
1738 
1757 template <int dim, int spacedim>
1758 class TriaAccessor<0, dim, spacedim>
1759 {
1760 public:
1766  static const unsigned int space_dimension = spacedim;
1767 
1773  static const unsigned int dimension = dim;
1774 
1780  static const unsigned int structure_dimension = 0;
1781 
1785  using AccessorData = void;
1786 
1792  const unsigned int vertex_index);
1793 
1799  TriaAccessor(const Triangulation<dim, spacedim> *tria = nullptr,
1800  const int level = 0,
1801  const int index = 0,
1802  const AccessorData * = nullptr);
1803 
1807  template <int structdim2, int dim2, int spacedim2>
1809 
1813  template <int structdim2, int dim2, int spacedim2>
1815 
1820  state() const;
1821 
1826  static int
1827  level();
1828 
1833  int
1834  index() const;
1835 
1841  get_triangulation() const;
1842 
1852  void
1853  operator++();
1854 
1858  void
1859  operator--();
1863  bool
1864  operator==(const TriaAccessor &) const;
1865 
1869  bool
1870  operator!=(const TriaAccessor &) const;
1871 
1899  unsigned int
1900  vertex_index(const unsigned int i = 0) const;
1901 
1907  Point<spacedim> &
1908  vertex(const unsigned int i = 0) const;
1909 
1914  typename ::internal::TriangulationImplementation::
1915  Iterators<dim, spacedim>::line_iterator static line(const unsigned int);
1916 
1920  static unsigned int
1921  line_index(const unsigned int i);
1922 
1926  static typename ::internal::TriangulationImplementation::
1927  Iterators<dim, spacedim>::quad_iterator
1928  quad(const unsigned int i);
1929 
1933  static unsigned int
1934  quad_index(const unsigned int i);
1935 
1951  double
1952  diameter() const;
1953 
1961  double
1962  extent_in_direction(const unsigned int axis) const;
1963 
1972  center(const bool respect_manifold = false,
1973  const bool interpolate_from_surrounding = false) const;
1974 
1983  double
1984  measure() const;
1999  static bool
2000  face_orientation(const unsigned int face);
2001 
2005  static bool
2006  face_flip(const unsigned int face);
2007 
2011  static bool
2012  face_rotation(const unsigned int face);
2013 
2017  static bool
2018  line_orientation(const unsigned int line);
2019 
2034  static bool
2035  has_children();
2036 
2041  static unsigned int
2042  n_children();
2043 
2048  static unsigned int
2049  number_of_children();
2050 
2054  static unsigned int
2055  max_refinement_depth();
2056 
2060  static unsigned int
2061  child_iterator_to_index(const TriaIterator<TriaAccessor<0, dim, spacedim>> &);
2062 
2067  child(const unsigned int);
2068 
2073  isotropic_child(const unsigned int);
2074 
2078  static RefinementCase<0>
2079  refinement_case();
2080 
2084  static int
2085  child_index(const unsigned int i);
2086 
2090  static int
2091  isotropic_child_index(const unsigned int i);
2099  bool
2100  used() const;
2101 
2102 protected:
2110  void
2111  copy_from(const TriaAccessor &);
2112 
2121  bool
2122  operator<(const TriaAccessor &other) const;
2123 
2128 
2132  unsigned int global_vertex_index;
2133 
2134 private:
2135  template <typename Accessor>
2136  friend class TriaRawIterator;
2137  template <typename Accessor>
2138  friend class TriaIterator;
2139  template <typename Accessor>
2140  friend class TriaActiveIterator;
2141 };
2142 
2143 
2144 
2161 template <int spacedim>
2162 class TriaAccessor<0, 1, spacedim>
2163 {
2164 public:
2170  static const unsigned int space_dimension = spacedim;
2171 
2177  static const unsigned int dimension = 1;
2178 
2184  static const unsigned int structure_dimension = 0;
2185 
2189  using AccessorData = void;
2190 
2196  {
2208  right_vertex
2209  };
2210 
2223  const VertexKind vertex_kind,
2224  const unsigned int vertex_index);
2225 
2231  TriaAccessor(const Triangulation<1, spacedim> *tria = nullptr,
2232  const int = 0,
2233  const int = 0,
2234  const AccessorData * = nullptr);
2235 
2239  template <int structdim2, int dim2, int spacedim2>
2241 
2245  template <int structdim2, int dim2, int spacedim2>
2247 
2252  void
2253  copy_from(const TriaAccessor &);
2254 
2261  state();
2262 
2267  static int
2268  level();
2269 
2274  int
2275  index() const;
2276 
2282  get_triangulation() const;
2283 
2294  void
2295  operator++() const;
2296 
2301  void
2302  operator--() const;
2306  bool
2307  operator==(const TriaAccessor &) const;
2308 
2312  bool
2313  operator!=(const TriaAccessor &) const;
2314 
2323  bool
2324  operator<(const TriaAccessor &other) const;
2325 
2352  unsigned int
2353  vertex_index(const unsigned int i = 0) const;
2354 
2360  Point<spacedim> &
2361  vertex(const unsigned int i = 0) const;
2362 
2368  center() const;
2369 
2374  typename ::internal::TriangulationImplementation::
2375  Iterators<1, spacedim>::line_iterator static line(const unsigned int);
2376 
2383  static unsigned int
2384  line_index(const unsigned int i);
2385 
2389  static typename ::internal::TriangulationImplementation::
2390  Iterators<1, spacedim>::quad_iterator
2391  quad(const unsigned int i);
2392 
2399  static unsigned int
2400  quad_index(const unsigned int i);
2401 
2411  bool
2412  at_boundary() const;
2413 
2429  boundary_id() const;
2430 
2434  const Manifold<1, spacedim> &
2435  get_manifold() const;
2436 
2444  manifold_id() const;
2445 
2446 
2457  static bool
2458  face_orientation(const unsigned int face);
2459 
2463  static bool
2464  face_flip(const unsigned int face);
2465 
2469  static bool
2470  face_rotation(const unsigned int face);
2471 
2475  static bool
2476  line_orientation(const unsigned int line);
2477 
2492  static bool
2493  has_children();
2494 
2499  static unsigned int
2500  n_children();
2501 
2506  static unsigned int
2507  number_of_children();
2508 
2512  static unsigned int
2513  max_refinement_depth();
2514 
2518  static unsigned int
2519  child_iterator_to_index(const TriaIterator<TriaAccessor<0, 1, spacedim>> &);
2520 
2525  child(const unsigned int);
2526 
2531  isotropic_child(const unsigned int);
2532 
2536  static RefinementCase<0>
2537  refinement_case();
2538 
2542  static int
2543  child_index(const unsigned int i);
2544 
2548  static int
2549  isotropic_child_index(const unsigned int i);
2582  void
2583  set_boundary_id(const types::boundary_id);
2584 
2591  void
2592  set_manifold_id(const types::manifold_id);
2593 
2605  void
2606  set_all_boundary_ids(const types::boundary_id);
2607 
2619  void
2620  set_all_manifold_ids(const types::manifold_id);
2628  bool
2629  used() const;
2630 
2631 protected:
2636 
2642 
2646  unsigned int global_vertex_index;
2647 };
2648 
2649 
2650 
2666 template <int dim, int spacedim = dim>
2667 class CellAccessor : public TriaAccessor<dim, dim, spacedim>
2668 {
2669 public:
2674 
2679 
2690  CellAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
2691  const int level = -1,
2692  const int index = -1,
2693  const AccessorData * local_data = nullptr);
2694 
2698  CellAccessor(const TriaAccessor<dim, dim, spacedim> &cell_accessor);
2699 
2712  template <int structdim2, int dim2, int spacedim2>
2714 
2719  template <int structdim2, int dim2, int spacedim2>
2721 
2731  void
2732  operator=(const CellAccessor<dim, spacedim> &) = delete;
2733 
2750  child(const unsigned int i) const;
2751 
2755  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>
2756  face(const unsigned int i) const;
2757 
2762  unsigned int
2763  face_iterator_to_index(
2765 
2769  std::array<TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>,
2771  face_iterators() const;
2772 
2782  unsigned int
2783  face_index(const unsigned int i) const;
2784 
2833  neighbor_child_on_subface(const unsigned int face_no,
2834  const unsigned int subface_no) const;
2835 
2858  neighbor(const unsigned int i) const;
2859 
2864  int
2865  neighbor_index(const unsigned int i) const;
2866 
2871  int
2872  neighbor_level(const unsigned int i) const;
2873 
2885  unsigned int
2886  neighbor_of_neighbor(const unsigned int neighbor) const;
2887 
2898  bool
2899  neighbor_is_coarser(const unsigned int neighbor) const;
2900 
2915  std::pair<unsigned int, unsigned int>
2916  neighbor_of_coarser_neighbor(const unsigned int neighbor) const;
2917 
2924  unsigned int
2925  neighbor_face_no(const unsigned int neighbor) const;
2926 
2930  static bool
2931  is_level_cell();
2932 
2946  bool
2947  has_periodic_neighbor(const unsigned int i) const;
2948 
2966  periodic_neighbor(const unsigned int i) const;
2967 
2976  neighbor_or_periodic_neighbor(const unsigned int i) const;
2977 
2993  periodic_neighbor_child_on_subface(const unsigned int face_no,
2994  const unsigned int subface_no) const;
2995 
3006  std::pair<unsigned int, unsigned int>
3007  periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const;
3008 
3014  int
3015  periodic_neighbor_index(const unsigned int i) const;
3016 
3022  int
3023  periodic_neighbor_level(const unsigned int i) const;
3024 
3039  unsigned int
3040  periodic_neighbor_of_periodic_neighbor(const unsigned int i) const;
3041 
3047  unsigned int
3048  periodic_neighbor_face_no(const unsigned int i) const;
3049 
3056  bool
3057  periodic_neighbor_is_coarser(const unsigned int i) const;
3058 
3075  bool
3076  at_boundary(const unsigned int i) const;
3077 
3086  bool
3087  at_boundary() const;
3088 
3096  bool
3097  has_boundary_lines() const;
3124  refine_flag_set() const;
3125 
3143  void
3144  set_refine_flag(const RefinementCase<dim> ref_case =
3146 
3150  void
3151  clear_refine_flag() const;
3152 
3160  bool
3161  flag_for_face_refinement(
3162  const unsigned int face_no,
3163  const RefinementCase<dim - 1> &face_refinement_case =
3165 
3171  bool
3172  flag_for_line_refinement(const unsigned int line_no) const;
3173 
3183  subface_case(const unsigned int face_no) const;
3184 
3188  bool
3189  coarsen_flag_set() const;
3190 
3195  void
3196  set_coarsen_flag() const;
3197 
3201  void
3202  clear_coarsen_flag() const;
3226  material_id() const;
3227 
3239  void
3240  set_material_id(const types::material_id new_material_id) const;
3241 
3250  void
3251  recursively_set_material_id(const types::material_id new_material_id) const;
3278  subdomain_id() const;
3279 
3295  void
3296  set_subdomain_id(const types::subdomain_id new_subdomain_id) const;
3297 
3307  level_subdomain_id() const;
3308 
3313  void
3314  set_level_subdomain_id(
3315  const types::subdomain_id new_level_subdomain_id) const;
3316 
3317 
3333  void
3334  recursively_set_subdomain_id(
3335  const types::subdomain_id new_subdomain_id) const;
3353  bool
3354  direction_flag() const;
3355 
3381  unsigned int
3382  active_cell_index() const;
3383 
3391  int
3392  parent_index() const;
3393 
3400  parent() const;
3401 
3427  bool
3428  active() const;
3429 
3438  bool
3439  is_active() const;
3440 
3460  bool
3461  is_locally_owned() const;
3462 
3467  bool
3468  is_locally_owned_on_level() const;
3469 
3493  bool
3494  is_ghost() const;
3495 
3522  bool
3523  is_artificial() const;
3524 
3538  bool
3539  point_inside(const Point<spacedim> &p) const;
3540 
3549  void
3550  set_neighbor(const unsigned int i,
3551  const TriaIterator<CellAccessor<dim, spacedim>> &pointer) const;
3552 
3566  CellId
3567  id() const;
3568 
3577  DeclException0(ExcRefineCellNotActive);
3581  DeclException0(ExcCellFlaggedForRefinement);
3585  DeclException0(ExcCellFlaggedForCoarsening);
3586 
3587 protected:
3603  unsigned int
3604  neighbor_of_neighbor_internal(const unsigned int neighbor) const;
3605 
3611  template <int dim_, int spacedim_>
3612  bool
3613  point_inside_codim(const Point<spacedim_> &p) const;
3614 
3615 
3616 
3617 private:
3622  void
3623  set_active_cell_index(const unsigned int active_cell_index);
3624 
3628  void
3629  set_parent(const unsigned int parent_index);
3630 
3637  void
3638  set_direction_flag(const bool new_direction_flag) const;
3639 
3640  template <int, int>
3641  friend class Triangulation;
3642 
3643  friend struct ::internal::TriangulationImplementation::Implementation;
3644 };
3645 
3646 
3647 
3648 /* ----- declaration of explicit specializations and general templates ----- */
3649 
3650 
3651 template <int structdim, int dim, int spacedim>
3652 template <typename OtherAccessor>
3654  const OtherAccessor &)
3655 {
3656  Assert(false,
3657  ExcMessage("You are attempting an illegal conversion between "
3658  "iterator/accessor types. The constructor you call "
3659  "only exists to make certain template constructs "
3660  "easier to write as dimension independent code but "
3661  "the conversion is not valid in the current context."));
3662 }
3663 
3664 
3665 
3666 template <int structdim, int dim, int spacedim>
3667 template <int structdim2, int dim2, int spacedim2>
3670 {
3671  Assert(false,
3672  ExcMessage("You are attempting an illegal conversion between "
3673  "iterator/accessor types. The constructor you call "
3674  "only exists to make certain template constructs "
3675  "easier to write as dimension independent code but "
3676  "the conversion is not valid in the current context."));
3677 }
3678 
3679 
3680 
3681 template <int dim, int spacedim>
3682 template <int structdim2, int dim2, int spacedim2>
3685 {
3686  Assert(false,
3687  ExcMessage("You are attempting an illegal conversion between "
3688  "iterator/accessor types. The constructor you call "
3689  "only exists to make certain template constructs "
3690  "easier to write as dimension independent code but "
3691  "the conversion is not valid in the current context."));
3692 }
3693 
3694 
3695 
3696 template <int structdim, int dim, int spacedim>
3697 template <int structdim2, int dim2, int spacedim2>
3700 {
3701  Assert(false,
3702  ExcMessage("You are attempting an illegal conversion between "
3703  "iterator/accessor types. The constructor you call "
3704  "only exists to make certain template constructs "
3705  "easier to write as dimension independent code but "
3706  "the conversion is not valid in the current context."));
3707 }
3708 
3709 
3710 
3711 template <int dim, int spacedim>
3712 template <int structdim2, int dim2, int spacedim2>
3715 {
3716  Assert(false,
3717  ExcMessage("You are attempting an illegal conversion between "
3718  "iterator/accessor types. The constructor you call "
3719  "only exists to make certain template constructs "
3720  "easier to write as dimension independent code but "
3721  "the conversion is not valid in the current context."));
3722 }
3723 
3724 
3725 #ifndef DOXYGEN
3726 
3727 template <>
3728 bool
3730 template <>
3731 bool
3733 template <>
3734 bool
3736 template <>
3737 bool
3739 template <>
3740 bool
3742 template <>
3743 bool
3745 // -------------------------------------------------------------------
3746 
3747 template <>
3748 void
3750 
3751 #endif // DOXYGEN
3752 
3754 
3755 // include more templates in debug and optimized mode
3756 #include "tria_accessor.templates.h"
3757 
3758 #endif
unsigned int manifold_id
Definition: types.h:141
static ::ExceptionBase & ExcCellHasNoChildren()
InvalidAccessor(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
static ::ExceptionBase & ExcNeighborIsNotCoarser()
double diameter(const Triangulation< dim, spacedim > &tria)
Definition: grid_tools.cc:76
const Triangulation< dim, spacedim > * tria
static ::ExceptionBase & ExcNeighborIsCoarser()
static ::ExceptionBase & ExcDereferenceInvalidObject(AccessorType arg1)
bool operator!=(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
const Triangulation< 1, spacedim > * tria
unsigned int material_id
Definition: types.h:152
void set_all_manifold_ids(const types::manifold_id) const
bool operator<(const SynchronousIterators< Iterators > &a, const SynchronousIterators< Iterators > &b)
static ::ExceptionBase & ExcFacesHaveNoLevel()
bool operator==(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
static ::ExceptionBase & ExcCellNotUsed()
static ::ExceptionBase & ExcMessage(std::string arg1)
unsigned int subdomain_id
Definition: types.h:43
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:515
const Triangulation< dim, spacedim > * tria
static ::ExceptionBase & ExcSetOnlyEvenChildren(int arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1403
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:493
#define DeclException0(Exception0)
Definition: exceptions.h:470
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:361
unsigned int level
Definition: grid_out.cc:4355
static ::ExceptionBase & ExcCellNotActive()
SynchronousIterators< Iterators > operator++(SynchronousIterators< Iterators > &a)
static ::ExceptionBase & ExcNoPeriodicNeighbor()
static ::ExceptionBase & ExcCantSetChildren(int arg1)
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
Definition: utilities.cc:188
Definition: cell_id.h:69
CellAccessor(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
static ::ExceptionBase & ExcCellHasNoParent()
bool point_inside(const Point< spacedim > &p) const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:360
Point< 3 > center
TriaAccessor(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
Iterator reached end of container.
Iterator points to a valid object.
#define DEAL_II_DEPRECATED
Definition: config.h:150
unsigned int boundary_id
Definition: types.h:129
typename ::internal::TriaAccessorImplementation::PresentLevelType< structdim, dim >::type present_level
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcCantCompareIterators()