Reference documentation for deal.II version Git 4abc4a1666 2020-07-04 19:58:34 +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 
1611  unsigned int
1612  n_vertices() const;
1613 
1617  unsigned int
1618  n_lines() const;
1619 
1625  unsigned int
1626  n_faces() const;
1627 
1633  vertex_indices() const;
1634 
1640  line_indices() const;
1641 
1649  face_indices() const;
1650 
1656 private:
1661  void
1662  set_boundary_id_internal(const types::boundary_id id) const;
1663 
1671  void
1672  set_bounding_object_indices(
1673  const std::initializer_list<int> &new_indices) const;
1674 
1678  void
1679  set_bounding_object_indices(
1680  const std::initializer_list<unsigned int> &new_indices) const;
1681 
1689  void
1690  set_line_orientation(const unsigned int line, const bool orientation) const;
1691 
1702  void
1703  set_face_orientation(const unsigned int face, const bool orientation) const;
1704 
1711  void
1712  set_face_flip(const unsigned int face, const bool flip) const;
1713 
1720  void
1721  set_face_rotation(const unsigned int face, const bool rotation) const;
1722 
1726  void
1727  set_used_flag() const;
1728 
1732  void
1733  clear_used_flag() const;
1734 
1743  void
1744  set_refinement_case(const RefinementCase<structdim> &ref_case) const;
1745 
1753  void
1754  clear_refinement_case() const;
1755 
1762  void
1763  set_children(const unsigned int i, const int index) const;
1764 
1769  void
1770  clear_children() const;
1771 
1772 private:
1773  template <int, int>
1774  friend class Triangulation;
1775 
1776  friend struct ::internal::TriangulationImplementation::Implementation;
1777  friend struct ::internal::TriaAccessorImplementation::Implementation;
1778 };
1779 
1780 
1781 
1800 template <int dim, int spacedim>
1801 class TriaAccessor<0, dim, spacedim>
1802 {
1803 public:
1809  static const unsigned int space_dimension = spacedim;
1810 
1816  static const unsigned int dimension = dim;
1817 
1823  static const unsigned int structure_dimension = 0;
1824 
1828  using AccessorData = void;
1829 
1835  const unsigned int vertex_index);
1836 
1842  TriaAccessor(const Triangulation<dim, spacedim> *tria = nullptr,
1843  const int level = 0,
1844  const int index = 0,
1845  const AccessorData * = nullptr);
1846 
1850  template <int structdim2, int dim2, int spacedim2>
1852 
1856  template <int structdim2, int dim2, int spacedim2>
1858 
1863  state() const;
1864 
1869  static int
1870  level();
1871 
1876  int
1877  index() const;
1878 
1884  get_triangulation() const;
1885 
1895  void
1896  operator++();
1897 
1901  void
1902  operator--();
1906  bool
1907  operator==(const TriaAccessor &) const;
1908 
1912  bool
1913  operator!=(const TriaAccessor &) const;
1914 
1942  unsigned int
1943  vertex_index(const unsigned int i = 0) const;
1944 
1950  Point<spacedim> &
1951  vertex(const unsigned int i = 0) const;
1952 
1957  typename ::internal::TriangulationImplementation::
1958  Iterators<dim, spacedim>::line_iterator static line(const unsigned int);
1959 
1963  static unsigned int
1964  line_index(const unsigned int i);
1965 
1969  static typename ::internal::TriangulationImplementation::
1970  Iterators<dim, spacedim>::quad_iterator
1971  quad(const unsigned int i);
1972 
1976  static unsigned int
1977  quad_index(const unsigned int i);
1978 
1994  double
1995  diameter() const;
1996 
2004  double
2005  extent_in_direction(const unsigned int axis) const;
2006 
2015  center(const bool respect_manifold = false,
2016  const bool interpolate_from_surrounding = false) const;
2017 
2026  double
2027  measure() const;
2042  static bool
2043  face_orientation(const unsigned int face);
2044 
2048  static bool
2049  face_flip(const unsigned int face);
2050 
2054  static bool
2055  face_rotation(const unsigned int face);
2056 
2060  static bool
2061  line_orientation(const unsigned int line);
2062 
2077  static bool
2078  has_children();
2079 
2084  static unsigned int
2085  n_children();
2086 
2091  static unsigned int
2092  number_of_children();
2093 
2097  static unsigned int
2098  max_refinement_depth();
2099 
2103  static unsigned int
2104  child_iterator_to_index(const TriaIterator<TriaAccessor<0, dim, spacedim>> &);
2105 
2110  child(const unsigned int);
2111 
2116  isotropic_child(const unsigned int);
2117 
2121  static RefinementCase<0>
2122  refinement_case();
2123 
2127  static int
2128  child_index(const unsigned int i);
2129 
2133  static int
2134  isotropic_child_index(const unsigned int i);
2142  bool
2143  used() const;
2144 
2145 protected:
2153  void
2154  copy_from(const TriaAccessor &);
2155 
2164  bool
2165  operator<(const TriaAccessor &other) const;
2166 
2171 
2175  unsigned int global_vertex_index;
2176 
2177 private:
2178  template <typename Accessor>
2179  friend class TriaRawIterator;
2180  template <typename Accessor>
2181  friend class TriaIterator;
2182  template <typename Accessor>
2183  friend class TriaActiveIterator;
2184 };
2185 
2186 
2187 
2204 template <int spacedim>
2205 class TriaAccessor<0, 1, spacedim>
2206 {
2207 public:
2213  static const unsigned int space_dimension = spacedim;
2214 
2220  static const unsigned int dimension = 1;
2221 
2227  static const unsigned int structure_dimension = 0;
2228 
2232  using AccessorData = void;
2233 
2239  {
2251  right_vertex
2252  };
2253 
2266  const VertexKind vertex_kind,
2267  const unsigned int vertex_index);
2268 
2274  TriaAccessor(const Triangulation<1, spacedim> *tria = nullptr,
2275  const int = 0,
2276  const int = 0,
2277  const AccessorData * = nullptr);
2278 
2282  template <int structdim2, int dim2, int spacedim2>
2284 
2288  template <int structdim2, int dim2, int spacedim2>
2290 
2295  void
2296  copy_from(const TriaAccessor &);
2297 
2304  state();
2305 
2310  static int
2311  level();
2312 
2317  int
2318  index() const;
2319 
2325  get_triangulation() const;
2326 
2337  void
2338  operator++() const;
2339 
2344  void
2345  operator--() const;
2349  bool
2350  operator==(const TriaAccessor &) const;
2351 
2355  bool
2356  operator!=(const TriaAccessor &) const;
2357 
2366  bool
2367  operator<(const TriaAccessor &other) const;
2368 
2395  unsigned int
2396  vertex_index(const unsigned int i = 0) const;
2397 
2403  Point<spacedim> &
2404  vertex(const unsigned int i = 0) const;
2405 
2411  center() const;
2412 
2417  typename ::internal::TriangulationImplementation::
2418  Iterators<1, spacedim>::line_iterator static line(const unsigned int);
2419 
2426  static unsigned int
2427  line_index(const unsigned int i);
2428 
2432  static typename ::internal::TriangulationImplementation::
2433  Iterators<1, spacedim>::quad_iterator
2434  quad(const unsigned int i);
2435 
2442  static unsigned int
2443  quad_index(const unsigned int i);
2444 
2454  bool
2455  at_boundary() const;
2456 
2472  boundary_id() const;
2473 
2477  const Manifold<1, spacedim> &
2478  get_manifold() const;
2479 
2487  manifold_id() const;
2488 
2489 
2500  static bool
2501  face_orientation(const unsigned int face);
2502 
2506  static bool
2507  face_flip(const unsigned int face);
2508 
2512  static bool
2513  face_rotation(const unsigned int face);
2514 
2518  static bool
2519  line_orientation(const unsigned int line);
2520 
2535  static bool
2536  has_children();
2537 
2542  static unsigned int
2543  n_children();
2544 
2549  static unsigned int
2550  number_of_children();
2551 
2555  static unsigned int
2556  max_refinement_depth();
2557 
2561  static unsigned int
2562  child_iterator_to_index(const TriaIterator<TriaAccessor<0, 1, spacedim>> &);
2563 
2568  child(const unsigned int);
2569 
2574  isotropic_child(const unsigned int);
2575 
2579  static RefinementCase<0>
2580  refinement_case();
2581 
2585  static int
2586  child_index(const unsigned int i);
2587 
2591  static int
2592  isotropic_child_index(const unsigned int i);
2625  void
2626  set_boundary_id(const types::boundary_id);
2627 
2634  void
2635  set_manifold_id(const types::manifold_id);
2636 
2648  void
2649  set_all_boundary_ids(const types::boundary_id);
2650 
2662  void
2663  set_all_manifold_ids(const types::manifold_id);
2671  bool
2672  used() const;
2673 
2677  unsigned int
2678  n_vertices() const;
2679 
2683  unsigned int
2684  n_lines() const;
2685 
2691  vertex_indices() const;
2692 
2698  line_indices() const;
2699 
2700 protected:
2705 
2711 
2715  unsigned int global_vertex_index;
2716 };
2717 
2718 
2719 
2735 template <int dim, int spacedim = dim>
2736 class CellAccessor : public TriaAccessor<dim, dim, spacedim>
2737 {
2738 public:
2743 
2748 
2759  CellAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
2760  const int level = -1,
2761  const int index = -1,
2762  const AccessorData * local_data = nullptr);
2763 
2767  CellAccessor(const TriaAccessor<dim, dim, spacedim> &cell_accessor);
2768 
2781  template <int structdim2, int dim2, int spacedim2>
2783 
2788  template <int structdim2, int dim2, int spacedim2>
2790 
2800  void
2801  operator=(const CellAccessor<dim, spacedim> &) = delete;
2802 
2819  child(const unsigned int i) const;
2820 
2824  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>
2825  face(const unsigned int i) const;
2826 
2831  unsigned int
2832  face_iterator_to_index(
2834 
2838  boost::container::small_vector<
2839  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>,
2841  face_iterators() const;
2842 
2852  unsigned int
2853  face_index(const unsigned int i) const;
2854 
2903  neighbor_child_on_subface(const unsigned int face_no,
2904  const unsigned int subface_no) const;
2905 
2928  neighbor(const unsigned int i) const;
2929 
2934  int
2935  neighbor_index(const unsigned int i) const;
2936 
2941  int
2942  neighbor_level(const unsigned int i) const;
2943 
2955  unsigned int
2956  neighbor_of_neighbor(const unsigned int neighbor) const;
2957 
2968  bool
2969  neighbor_is_coarser(const unsigned int neighbor) const;
2970 
2985  std::pair<unsigned int, unsigned int>
2986  neighbor_of_coarser_neighbor(const unsigned int neighbor) const;
2987 
2994  unsigned int
2995  neighbor_face_no(const unsigned int neighbor) const;
2996 
3000  static bool
3001  is_level_cell();
3002 
3016  bool
3017  has_periodic_neighbor(const unsigned int i) const;
3018 
3036  periodic_neighbor(const unsigned int i) const;
3037 
3046  neighbor_or_periodic_neighbor(const unsigned int i) const;
3047 
3063  periodic_neighbor_child_on_subface(const unsigned int face_no,
3064  const unsigned int subface_no) const;
3065 
3076  std::pair<unsigned int, unsigned int>
3077  periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const;
3078 
3084  int
3085  periodic_neighbor_index(const unsigned int i) const;
3086 
3092  int
3093  periodic_neighbor_level(const unsigned int i) const;
3094 
3109  unsigned int
3110  periodic_neighbor_of_periodic_neighbor(const unsigned int i) const;
3111 
3117  unsigned int
3118  periodic_neighbor_face_no(const unsigned int i) const;
3119 
3126  bool
3127  periodic_neighbor_is_coarser(const unsigned int i) const;
3128 
3145  bool
3146  at_boundary(const unsigned int i) const;
3147 
3156  bool
3157  at_boundary() const;
3158 
3166  bool
3167  has_boundary_lines() const;
3194  refine_flag_set() const;
3195 
3213  void
3214  set_refine_flag(const RefinementCase<dim> ref_case =
3216 
3220  void
3221  clear_refine_flag() const;
3222 
3230  bool
3231  flag_for_face_refinement(
3232  const unsigned int face_no,
3233  const RefinementCase<dim - 1> &face_refinement_case =
3235 
3241  bool
3242  flag_for_line_refinement(const unsigned int line_no) const;
3243 
3253  subface_case(const unsigned int face_no) const;
3254 
3258  bool
3259  coarsen_flag_set() const;
3260 
3265  void
3266  set_coarsen_flag() const;
3267 
3271  void
3272  clear_coarsen_flag() const;
3296  material_id() const;
3297 
3309  void
3310  set_material_id(const types::material_id new_material_id) const;
3311 
3320  void
3321  recursively_set_material_id(const types::material_id new_material_id) const;
3348  subdomain_id() const;
3349 
3365  void
3366  set_subdomain_id(const types::subdomain_id new_subdomain_id) const;
3367 
3377  level_subdomain_id() const;
3378 
3383  void
3384  set_level_subdomain_id(
3385  const types::subdomain_id new_level_subdomain_id) const;
3386 
3387 
3403  void
3404  recursively_set_subdomain_id(
3405  const types::subdomain_id new_subdomain_id) const;
3423  bool
3424  direction_flag() const;
3425 
3451  unsigned int
3452  active_cell_index() const;
3453 
3461  int
3462  parent_index() const;
3463 
3470  parent() const;
3471 
3497  bool
3498  active() const;
3499 
3508  bool
3509  is_active() const;
3510 
3530  bool
3531  is_locally_owned() const;
3532 
3537  bool
3538  is_locally_owned_on_level() const;
3539 
3563  bool
3564  is_ghost() const;
3565 
3592  bool
3593  is_artificial() const;
3594 
3608  bool
3609  point_inside(const Point<spacedim> &p) const;
3610 
3619  void
3620  set_neighbor(const unsigned int i,
3621  const TriaIterator<CellAccessor<dim, spacedim>> &pointer) const;
3622 
3636  CellId
3637  id() const;
3638 
3647  DeclException0(ExcRefineCellNotActive);
3651  DeclException0(ExcCellFlaggedForRefinement);
3655  DeclException0(ExcCellFlaggedForCoarsening);
3656 
3657 protected:
3673  unsigned int
3674  neighbor_of_neighbor_internal(const unsigned int neighbor) const;
3675 
3681  template <int dim_, int spacedim_>
3682  bool
3683  point_inside_codim(const Point<spacedim_> &p) const;
3684 
3685 
3686 
3687 private:
3692  void
3693  set_active_cell_index(const unsigned int active_cell_index);
3694 
3698  void
3699  set_parent(const unsigned int parent_index);
3700 
3707  void
3708  set_direction_flag(const bool new_direction_flag) const;
3709 
3710  template <int, int>
3711  friend class Triangulation;
3712 
3713  friend struct ::internal::TriangulationImplementation::Implementation;
3714 };
3715 
3716 
3717 
3718 /* ----- declaration of explicit specializations and general templates ----- */
3719 
3720 
3721 template <int structdim, int dim, int spacedim>
3722 template <typename OtherAccessor>
3724  const OtherAccessor &)
3725 {
3726  Assert(false,
3727  ExcMessage("You are attempting an illegal conversion between "
3728  "iterator/accessor types. The constructor you call "
3729  "only exists to make certain template constructs "
3730  "easier to write as dimension independent code but "
3731  "the conversion is not valid in the current context."));
3732 }
3733 
3734 
3735 
3736 template <int structdim, int dim, int spacedim>
3737 template <int structdim2, int dim2, int spacedim2>
3740 {
3741  Assert(false,
3742  ExcMessage("You are attempting an illegal conversion between "
3743  "iterator/accessor types. The constructor you call "
3744  "only exists to make certain template constructs "
3745  "easier to write as dimension independent code but "
3746  "the conversion is not valid in the current context."));
3747 }
3748 
3749 
3750 
3751 template <int dim, int spacedim>
3752 template <int structdim2, int dim2, int spacedim2>
3755 {
3756  Assert(false,
3757  ExcMessage("You are attempting an illegal conversion between "
3758  "iterator/accessor types. The constructor you call "
3759  "only exists to make certain template constructs "
3760  "easier to write as dimension independent code but "
3761  "the conversion is not valid in the current context."));
3762 }
3763 
3764 
3765 
3766 template <int structdim, int dim, int spacedim>
3767 template <int structdim2, int dim2, int spacedim2>
3770 {
3771  Assert(false,
3772  ExcMessage("You are attempting an illegal conversion between "
3773  "iterator/accessor types. The constructor you call "
3774  "only exists to make certain template constructs "
3775  "easier to write as dimension independent code but "
3776  "the conversion is not valid in the current context."));
3777 }
3778 
3779 
3780 
3781 template <int dim, int spacedim>
3782 template <int structdim2, int dim2, int spacedim2>
3785 {
3786  Assert(false,
3787  ExcMessage("You are attempting an illegal conversion between "
3788  "iterator/accessor types. The constructor you call "
3789  "only exists to make certain template constructs "
3790  "easier to write as dimension independent code but "
3791  "the conversion is not valid in the current context."));
3792 }
3793 
3794 
3795 #ifndef DOXYGEN
3796 
3797 template <>
3798 bool
3800 template <>
3801 bool
3803 template <>
3804 bool
3806 template <>
3807 bool
3809 template <>
3810 bool
3812 template <>
3813 bool
3815 // -------------------------------------------------------------------
3816 
3817 template <>
3818 void
3820 
3821 #endif // DOXYGEN
3822 
3824 
3825 // include more templates in debug and optimized mode
3826 #include "tria_accessor.templates.h"
3827 
3828 #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:78
const Triangulation< dim, spacedim > * tria
static ::ExceptionBase & ExcNeighborIsCoarser()
static ::ExceptionBase & ExcDereferenceInvalidObject(AccessorType arg1)
std::vector< unsigned int > vertex_indices
Definition: tria.cc:2244
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
boost::integer_range< IncrementableType > iota_view
Definition: iota_view.h:46
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:493
#define DeclException0(Exception0)
Definition: exceptions.h:470
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
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:362
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:152
unsigned int boundary_id
Definition: types.h:129
typename ::internal::TriaAccessorImplementation::PresentLevelType< structdim, dim >::type present_level
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcCantCompareIterators()