Reference documentation for deal.II version Git 87eb8ae23d 2021-01-27 12:59:39 +0100
\(\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>
31 
32 #include <boost/container/small_vector.hpp>
33 
34 #include <utility>
35 
36 
38 
39 // Forward declarations
40 #ifndef DOXYGEN
41 template <int dim, int spacedim>
42 class Triangulation;
43 template <typename Accessor>
44 class TriaRawIterator;
45 template <typename Accessor>
46 class TriaIterator;
47 template <typename Accessor>
48 class TriaActiveIterator;
49 
50 namespace parallel
51 {
52  template <int dim, int spacedim>
53  class TriangulationBase;
54 }
55 
56 template <int dim, int spacedim>
57 class Manifold;
58 
59 template <int dim, int spacedim>
60 class Mapping;
61 #endif
62 
63 namespace internal
64 {
65  namespace TriangulationImplementation
66  {
67  class TriaObjects;
68  struct Implementation;
69  struct ImplementationMixedMesh;
70  } // namespace TriangulationImplementation
71 
72  namespace TriaAccessorImplementation
73  {
74  struct Implementation;
75 
81  template <int structdim, int dim>
83  {
84  struct type
85  {
89  type() = default;
90 
94  type(const int level)
95  {
96  Assert(level == 0, ExcInternalError());
97  (void)level; // removes -Wunused-parameter warning in optimized mode
98  }
99 
103  operator int() const
104  {
105  return 0;
106  }
107 
108  void
109  operator++() const
110  {
111  Assert(false, ExcInternalError());
112  }
113 
114  void
115  operator--() const
116  {
117  Assert(false, ExcInternalError());
118  }
119  };
120  };
121 
122 
128  template <int dim>
129  struct PresentLevelType<dim, dim>
130  {
131  using type = int;
132  };
133  } // namespace TriaAccessorImplementation
134 } // namespace internal
135 template <int structdim, int dim, int spacedim>
137 template <int dim, int spacedim>
138 class TriaAccessor<0, dim, spacedim>;
139 template <int spacedim>
140 class TriaAccessor<0, 1, spacedim>;
141 
142 // note: the file tria_accessor.templates.h is included at the end of
143 // this file. this includes a lot of templates. originally, this was
144 // only done in debug mode, but led to cyclic reduction problems and
145 // so is now on by default.
146 
147 
152 {
157  "The operation you are attempting can only be performed for "
158  "(cell, face, or edge) iterators that point to valid "
159  "objects. These objects need not necessarily be active, "
160  "i.e., have no children, but they need to be part of a "
161  "triangulation. (The objects pointed to by an iterator "
162  "may -- after coarsening -- also be objects that used "
163  "to be part of a triangulation, but are now no longer "
164  "used. Their memory location may have been retained "
165  "for re-use upon the next mesh refinement, but is "
166  "currently unused.)");
177  "The operation you are attempting can only be performed for "
178  "(cell, face, or edge) iterators that point to 'active' "
179  "objects. 'Active' objects are those that do not have "
180  "children (in the case of cells), or that are part of "
181  "an active cell (in the case of faces or edges). However, "
182  "the object on which you are trying the current "
183  "operation is not 'active' in this sense.");
190  "The operation you are attempting can only be performed for "
191  "(cell, face, or edge) iterators that have children, "
192  "but the object on which you are trying the current "
193  "operation does not have any.");
201  "The operation you are attempting can only be performed for "
202  "(cell, face, or edge) iterators that have a parent object, "
203  "but the object on which you are trying the current "
204  "operation does not have one -- i.e., it is on the "
205  "coarsest level of the triangulation.");
210  int,
211  << "You can only set the child index if the cell does not "
212  << "currently have children registered; or you can clear it. "
213  << "The given index was " << arg1
214  << " (-1 means: clear children).");
218  template <typename AccessorType>
220  AccessorType,
221  << "You tried to dereference an iterator for which this "
222  << "is not possible. More information on this iterator: "
223  << "index=" << arg1.index() << ", state="
224  << (arg1.state() == IteratorState::valid ?
225  "valid" :
226  (arg1.state() == IteratorState::past_the_end ?
227  "past_the_end" :
228  "invalid")));
233  "Iterators can only be compared if they point to the same "
234  "triangulation, or if neither of them are associated "
235  "with a triangulation.");
236  // TODO: Write documentation!
241  // TODO: Write documentation!
261  // TODO: Write documentation!
267  int,
268  << "You can only set the child index of an even numbered child."
269  << "The number of the child given was " << arg1 << ".");
270 } // namespace TriaAccessorExceptions
271 
272 
298 template <int structdim, int dim, int spacedim = dim>
300 {
301 public:
307  static const unsigned int space_dimension = spacedim;
308 
314  static const unsigned int dimension = dim;
315 
321  static const unsigned int structure_dimension = structdim;
322 
332  void
333  operator=(const TriaAccessorBase *) = delete;
334 
335 protected:
341  using AccessorData = void;
342 
346  TriaAccessorBase(const Triangulation<dim, spacedim> *parent = nullptr,
347  const int level = -1,
348  const int index = -1,
349  const AccessorData * = nullptr);
350 
355 
363  void
364  copy_from(const TriaAccessorBase &);
365 
370  operator=(const TriaAccessorBase &);
371 
382  bool
383  operator<(const TriaAccessorBase &other) const;
384 
385 protected:
389  bool
390  operator==(const TriaAccessorBase &) const;
391 
395  bool
396  operator!=(const TriaAccessorBase &) const;
397 
411  void
412  operator++();
413 
421  void
422  operator--();
431  objects() const;
432 
433 public:
439  using LocalData = void *;
440 
464  int
465  level() const;
466 
493  int
494  index() const;
495 
501  state() const;
502 
508  get_triangulation() const;
509 
513 protected:
518  typename ::internal::TriaAccessorImplementation::
519  PresentLevelType<structdim, dim>::type present_level;
520 
526 
531 
532 private:
533  template <typename Accessor>
534  friend class TriaRawIterator;
535  template <typename Accessor>
536  friend class TriaIterator;
537  template <typename Accessor>
538  friend class TriaActiveIterator;
539 };
540 
541 
542 
563 template <int structdim, int dim, int spacedim = dim>
564 class InvalidAccessor : public TriaAccessorBase<structdim, dim, spacedim>
565 {
566 public:
570  using AccessorData =
572 
580  InvalidAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
581  const int level = -1,
582  const int index = -1,
583  const AccessorData * local_data = nullptr);
584 
593 
598  template <typename OtherAccessor>
599  InvalidAccessor(const OtherAccessor &);
600 
604  void
605  copy_from(const InvalidAccessor &);
606 
610  bool
611  operator==(const InvalidAccessor &) const;
612  bool
613  operator!=(const InvalidAccessor &) const;
614 
618  void
619  operator++() const;
620  void
621  operator--() const;
622 
627  bool
628  used() const;
629 
634  bool
635  has_children() const;
636 
641  manifold_id() const;
642 
646  unsigned int
647  user_index() const;
648 
652  void
653  set_user_index(const unsigned int p) const;
654 
658  void
659  set_manifold_id(const types::manifold_id) const;
660 
665  vertex(const unsigned int i) const;
666 
671  typename ::internal::TriangulationImplementation::
672  Iterators<dim, spacedim>::line_iterator
673  line(const unsigned int i) const;
674 
679  typename ::internal::TriangulationImplementation::
680  Iterators<dim, spacedim>::quad_iterator
681  quad(const unsigned int i) const;
682 };
683 
684 
685 
703 template <int structdim, int dim, int spacedim>
704 class TriaAccessor : public TriaAccessorBase<structdim, dim, spacedim>
705 {
706 public:
710  using AccessorData =
712 
716  TriaAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
717  const int level = -1,
718  const int index = -1,
719  const AccessorData * local_data = nullptr);
720 
733  template <int structdim2, int dim2, int spacedim2>
735 
740  template <int structdim2, int dim2, int spacedim2>
742 
752  void
753  operator=(const TriaAccessor &) = delete;
754 
761  bool
762  used() const;
763 
776  vertex_iterator(const unsigned int i) const;
777 
793  unsigned int
794  vertex_index(const unsigned int i) const;
795 
834  vertex(const unsigned int i) const;
835 
839  typename ::internal::TriangulationImplementation::
840  Iterators<dim, spacedim>::line_iterator
841  line(const unsigned int i) const;
842 
849  unsigned int
850  line_index(const unsigned int i) const;
851 
855  typename ::internal::TriangulationImplementation::
856  Iterators<dim, spacedim>::quad_iterator
857  quad(const unsigned int i) const;
858 
865  unsigned int
866  quad_index(const unsigned int i) const;
889  bool
890  face_orientation(const unsigned int face) const;
891 
901  bool
902  face_flip(const unsigned int face) const;
903 
913  bool
914  face_rotation(const unsigned int face) const;
915 
926  bool
927  line_orientation(const unsigned int line) const;
942  bool
943  has_children() const;
944 
949  unsigned int
950  n_children() const;
951 
965  unsigned int
966  number_of_children() const;
967 
981  unsigned int
982  max_refinement_depth() const;
983 
988  child(const unsigned int i) const;
989 
994  unsigned int
995  child_iterator_to_index(
997 
1007  isotropic_child(const unsigned int i) const;
1008 
1013  refinement_case() const;
1014 
1020  int
1021  child_index(const unsigned int i) const;
1022 
1028  int
1029  isotropic_child_index(const unsigned int i) const;
1052  boundary_id() const;
1053 
1083  void
1084  set_boundary_id(const types::boundary_id) const;
1085 
1116  void
1117  set_all_boundary_ids(const types::boundary_id) const;
1118 
1126  bool
1127  at_boundary() const;
1128 
1138  const Manifold<dim, spacedim> &
1139  get_manifold() const;
1140 
1162  manifold_id() const;
1163 
1181  void
1182  set_manifold_id(const types::manifold_id) const;
1183 
1197  void
1198  set_all_manifold_ids(const types::manifold_id) const;
1199 
1216  bool
1217  user_flag_set() const;
1218 
1224  void
1225  set_user_flag() const;
1226 
1232  void
1233  clear_user_flag() const;
1234 
1240  void
1241  recursively_set_user_flag() const;
1242 
1248  void
1249  recursively_clear_user_flag() const;
1250 
1256  void
1257  clear_user_data() const;
1258 
1270  void
1271  set_user_pointer(void *p) const;
1272 
1278  void
1279  clear_user_pointer() const;
1280 
1296  void *
1297  user_pointer() const;
1298 
1320  void
1321  recursively_set_user_pointer(void *p) const;
1322 
1329  void
1330  recursively_clear_user_pointer() const;
1331 
1341  void
1342  set_user_index(const unsigned int p) const;
1343 
1349  void
1350  clear_user_index() const;
1351 
1363  unsigned int
1364  user_index() const;
1365 
1383  void
1384  recursively_set_user_index(const unsigned int p) const;
1385 
1394  void
1395  recursively_clear_user_index() const;
1430  double
1431  diameter() const;
1432 
1459  std::pair<Point<spacedim>, double>
1460  enclosing_ball() const;
1461 
1471  bounding_box() const;
1472 
1482  double
1483  extent_in_direction(const unsigned int axis) const;
1484 
1488  double
1489  minimum_vertex_distance() const;
1490 
1505  intermediate_point(const Point<structdim> &coordinates) const;
1506 
1530  real_to_unit_cell_affine_approximation(const Point<spacedim> &point) const;
1531 
1567  center(const bool respect_manifold = false,
1568  const bool interpolate_from_surrounding = false) const;
1569 
1588  barycenter() const;
1589 
1615  double
1616  measure() const;
1617 
1632  bool
1633  is_translation_of(
1635 
1640  reference_cell_type() const;
1641 
1645  unsigned int
1646  n_vertices() const;
1647 
1651  unsigned int
1652  n_lines() const;
1653 
1659  unsigned int
1660  n_faces() const;
1661 
1667  vertex_indices() const;
1668 
1674  line_indices() const;
1675 
1683  face_indices() const;
1684 
1689 protected:
1695  reference_cell_info() const;
1696 
1697 private:
1702  void
1703  set_boundary_id_internal(const types::boundary_id id) const;
1704 
1712  void
1713  set_bounding_object_indices(
1714  const std::initializer_list<int> &new_indices) const;
1715 
1719  void
1720  set_bounding_object_indices(
1721  const std::initializer_list<unsigned int> &new_indices) const;
1722 
1730  void
1731  set_line_orientation(const unsigned int line, const bool orientation) const;
1732 
1743  void
1744  set_face_orientation(const unsigned int face, const bool orientation) const;
1745 
1752  void
1753  set_face_flip(const unsigned int face, const bool flip) const;
1754 
1761  void
1762  set_face_rotation(const unsigned int face, const bool rotation) const;
1763 
1767  void
1768  set_used_flag() const;
1769 
1773  void
1774  clear_used_flag() const;
1775 
1784  void
1785  set_refinement_case(const RefinementCase<structdim> &ref_case) const;
1786 
1794  void
1795  clear_refinement_case() const;
1796 
1803  void
1804  set_children(const unsigned int i, const int index) const;
1805 
1810  void
1811  clear_children() const;
1812 
1813 private:
1814  template <int, int>
1815  friend class Triangulation;
1816 
1817  friend struct ::internal::TriangulationImplementation::Implementation;
1818  friend struct ::internal::TriangulationImplementation::
1819  ImplementationMixedMesh;
1820  friend struct ::internal::TriaAccessorImplementation::Implementation;
1821 };
1822 
1823 
1824 
1843 template <int dim, int spacedim>
1844 class TriaAccessor<0, dim, spacedim>
1845 {
1846 public:
1852  static const unsigned int space_dimension = spacedim;
1853 
1859  static const unsigned int dimension = dim;
1860 
1866  static const unsigned int structure_dimension = 0;
1867 
1871  using AccessorData = void;
1872 
1878  const unsigned int vertex_index);
1879 
1885  TriaAccessor(const Triangulation<dim, spacedim> *tria = nullptr,
1886  const int level = 0,
1887  const int index = 0,
1888  const AccessorData * = nullptr);
1889 
1893  template <int structdim2, int dim2, int spacedim2>
1895 
1899  template <int structdim2, int dim2, int spacedim2>
1901 
1906  state() const;
1907 
1912  static int
1913  level();
1914 
1919  int
1920  index() const;
1921 
1927  get_triangulation() const;
1928 
1938  void
1939  operator++();
1940 
1944  void
1945  operator--();
1949  bool
1950  operator==(const TriaAccessor &) const;
1951 
1955  bool
1956  operator!=(const TriaAccessor &) const;
1957 
1985  unsigned int
1986  vertex_index(const unsigned int i = 0) const;
1987 
1993  Point<spacedim> &
1994  vertex(const unsigned int i = 0) const;
1995 
2000  typename ::internal::TriangulationImplementation::
2001  Iterators<dim, spacedim>::line_iterator static line(const unsigned int);
2002 
2006  static unsigned int
2007  line_index(const unsigned int i);
2008 
2012  static typename ::internal::TriangulationImplementation::
2013  Iterators<dim, spacedim>::quad_iterator
2014  quad(const unsigned int i);
2015 
2019  static unsigned int
2020  quad_index(const unsigned int i);
2021 
2037  double
2038  diameter() const;
2039 
2047  double
2048  extent_in_direction(const unsigned int axis) const;
2049 
2058  center(const bool respect_manifold = false,
2059  const bool interpolate_from_surrounding = false) const;
2060 
2069  double
2070  measure() const;
2085  static bool
2086  face_orientation(const unsigned int face);
2087 
2091  static bool
2092  face_flip(const unsigned int face);
2093 
2097  static bool
2098  face_rotation(const unsigned int face);
2099 
2103  static bool
2104  line_orientation(const unsigned int line);
2105 
2120  static bool
2121  has_children();
2122 
2127  static unsigned int
2128  n_children();
2129 
2134  static unsigned int
2135  number_of_children();
2136 
2140  static unsigned int
2141  max_refinement_depth();
2142 
2146  static unsigned int
2147  child_iterator_to_index(const TriaIterator<TriaAccessor<0, dim, spacedim>> &);
2148 
2153  child(const unsigned int);
2154 
2159  isotropic_child(const unsigned int);
2160 
2164  static RefinementCase<0>
2165  refinement_case();
2166 
2170  static int
2171  child_index(const unsigned int i);
2172 
2176  static int
2177  isotropic_child_index(const unsigned int i);
2185  bool
2186  used() const;
2187 
2188 protected:
2196  void
2197  copy_from(const TriaAccessor &);
2198 
2207  bool
2208  operator<(const TriaAccessor &other) const;
2209 
2214 
2218  unsigned int global_vertex_index;
2219 
2220 private:
2221  template <typename Accessor>
2222  friend class TriaRawIterator;
2223  template <typename Accessor>
2224  friend class TriaIterator;
2225  template <typename Accessor>
2226  friend class TriaActiveIterator;
2227 };
2228 
2229 
2230 
2247 template <int spacedim>
2248 class TriaAccessor<0, 1, spacedim>
2249 {
2250 public:
2256  static const unsigned int space_dimension = spacedim;
2257 
2263  static const unsigned int dimension = 1;
2264 
2270  static const unsigned int structure_dimension = 0;
2271 
2275  using AccessorData = void;
2276 
2282  {
2294  right_vertex
2295  };
2296 
2309  const VertexKind vertex_kind,
2310  const unsigned int vertex_index);
2311 
2317  TriaAccessor(const Triangulation<1, spacedim> *tria = nullptr,
2318  const int = 0,
2319  const int = 0,
2320  const AccessorData * = nullptr);
2321 
2325  template <int structdim2, int dim2, int spacedim2>
2327 
2331  template <int structdim2, int dim2, int spacedim2>
2333 
2338  void
2339  copy_from(const TriaAccessor &);
2340 
2347  state();
2348 
2353  static int
2354  level();
2355 
2360  int
2361  index() const;
2362 
2368  get_triangulation() const;
2369 
2380  void
2381  operator++() const;
2382 
2387  void
2388  operator--() const;
2392  bool
2393  operator==(const TriaAccessor &) const;
2394 
2398  bool
2399  operator!=(const TriaAccessor &) const;
2400 
2409  bool
2410  operator<(const TriaAccessor &other) const;
2411 
2438  unsigned int
2439  vertex_index(const unsigned int i = 0) const;
2440 
2446  Point<spacedim> &
2447  vertex(const unsigned int i = 0) const;
2448 
2454  center() const;
2455 
2460  typename ::internal::TriangulationImplementation::
2461  Iterators<1, spacedim>::line_iterator static line(const unsigned int);
2462 
2469  static unsigned int
2470  line_index(const unsigned int i);
2471 
2475  static typename ::internal::TriangulationImplementation::
2476  Iterators<1, spacedim>::quad_iterator
2477  quad(const unsigned int i);
2478 
2485  static unsigned int
2486  quad_index(const unsigned int i);
2487 
2497  bool
2498  at_boundary() const;
2499 
2515  boundary_id() const;
2516 
2520  const Manifold<1, spacedim> &
2521  get_manifold() const;
2522 
2530  manifold_id() const;
2531 
2532 
2543  static bool
2544  face_orientation(const unsigned int face);
2545 
2549  static bool
2550  face_flip(const unsigned int face);
2551 
2555  static bool
2556  face_rotation(const unsigned int face);
2557 
2561  static bool
2562  line_orientation(const unsigned int line);
2563 
2578  static bool
2579  has_children();
2580 
2585  static unsigned int
2586  n_children();
2587 
2592  static unsigned int
2593  number_of_children();
2594 
2598  static unsigned int
2599  max_refinement_depth();
2600 
2604  static unsigned int
2605  child_iterator_to_index(const TriaIterator<TriaAccessor<0, 1, spacedim>> &);
2606 
2611  child(const unsigned int);
2612 
2617  isotropic_child(const unsigned int);
2618 
2622  static RefinementCase<0>
2623  refinement_case();
2624 
2628  static int
2629  child_index(const unsigned int i);
2630 
2634  static int
2635  isotropic_child_index(const unsigned int i);
2668  void
2669  set_boundary_id(const types::boundary_id);
2670 
2677  void
2678  set_manifold_id(const types::manifold_id);
2679 
2691  void
2692  set_all_boundary_ids(const types::boundary_id);
2693 
2705  void
2706  set_all_manifold_ids(const types::manifold_id);
2714  bool
2715  used() const;
2716 
2721  reference_cell_type() const;
2722 
2726  unsigned int
2727  n_vertices() const;
2728 
2732  unsigned int
2733  n_lines() const;
2734 
2740  vertex_indices() const;
2741 
2747  line_indices() const;
2748 
2749 protected:
2754 
2760 
2764  unsigned int global_vertex_index;
2765 };
2766 
2767 
2768 
2784 template <int dim, int spacedim = dim>
2785 class CellAccessor : public TriaAccessor<dim, dim, spacedim>
2786 {
2787 public:
2792 
2797 
2808  CellAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
2809  const int level = -1,
2810  const int index = -1,
2811  const AccessorData * local_data = nullptr);
2812 
2816  CellAccessor(const TriaAccessor<dim, dim, spacedim> &cell_accessor);
2817 
2830  template <int structdim2, int dim2, int spacedim2>
2832 
2837  template <int structdim2, int dim2, int spacedim2>
2839 
2849  void
2850  operator=(const CellAccessor<dim, spacedim> &) = delete;
2851 
2868  child(const unsigned int i) const;
2869 
2873  boost::container::small_vector<TriaIterator<CellAccessor<dim, spacedim>>,
2875  child_iterators() const;
2876 
2880  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>
2881  face(const unsigned int i) const;
2882 
2887  unsigned int
2888  face_iterator_to_index(
2890 
2894  boost::container::small_vector<
2895  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>,
2897  face_iterators() const;
2898 
2908  unsigned int
2909  face_index(const unsigned int i) const;
2910 
2959  neighbor_child_on_subface(const unsigned int face_no,
2960  const unsigned int subface_no) const;
2961 
2987  neighbor(const unsigned int face_no) const;
2988 
2996  int
2997  neighbor_index(const unsigned int face_no) const;
2998 
3006  int
3007  neighbor_level(const unsigned int face_no) const;
3008 
3020  unsigned int
3021  neighbor_of_neighbor(const unsigned int face_no) const;
3022 
3033  bool
3034  neighbor_is_coarser(const unsigned int face_no) const;
3035 
3050  std::pair<unsigned int, unsigned int>
3051  neighbor_of_coarser_neighbor(const unsigned int neighbor) const;
3052 
3059  unsigned int
3060  neighbor_face_no(const unsigned int neighbor) const;
3061 
3065  static bool
3066  is_level_cell();
3067 
3081  bool
3082  has_periodic_neighbor(const unsigned int i) const;
3083 
3101  periodic_neighbor(const unsigned int i) const;
3102 
3111  neighbor_or_periodic_neighbor(const unsigned int i) const;
3112 
3128  periodic_neighbor_child_on_subface(const unsigned int face_no,
3129  const unsigned int subface_no) const;
3130 
3141  std::pair<unsigned int, unsigned int>
3142  periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const;
3143 
3149  int
3150  periodic_neighbor_index(const unsigned int i) const;
3151 
3157  int
3158  periodic_neighbor_level(const unsigned int i) const;
3159 
3174  unsigned int
3175  periodic_neighbor_of_periodic_neighbor(const unsigned int i) const;
3176 
3182  unsigned int
3183  periodic_neighbor_face_no(const unsigned int i) const;
3184 
3191  bool
3192  periodic_neighbor_is_coarser(const unsigned int i) const;
3193 
3210  bool
3211  at_boundary(const unsigned int i) const;
3212 
3221  bool
3222  at_boundary() const;
3223 
3231  bool
3232  has_boundary_lines() const;
3259  refine_flag_set() const;
3260 
3278  void
3279  set_refine_flag(const RefinementCase<dim> ref_case =
3281 
3285  void
3286  clear_refine_flag() const;
3287 
3295  bool
3296  flag_for_face_refinement(
3297  const unsigned int face_no,
3298  const RefinementCase<dim - 1> &face_refinement_case =
3300 
3306  bool
3307  flag_for_line_refinement(const unsigned int line_no) const;
3308 
3318  subface_case(const unsigned int face_no) const;
3319 
3323  bool
3324  coarsen_flag_set() const;
3325 
3330  void
3331  set_coarsen_flag() const;
3332 
3336  void
3337  clear_coarsen_flag() const;
3361  material_id() const;
3362 
3374  void
3375  set_material_id(const types::material_id new_material_id) const;
3376 
3385  void
3386  recursively_set_material_id(const types::material_id new_material_id) const;
3413  subdomain_id() const;
3414 
3430  void
3431  set_subdomain_id(const types::subdomain_id new_subdomain_id) const;
3432 
3442  level_subdomain_id() const;
3443 
3448  void
3449  set_level_subdomain_id(
3450  const types::subdomain_id new_level_subdomain_id) const;
3451 
3452 
3468  void
3469  recursively_set_subdomain_id(
3470  const types::subdomain_id new_subdomain_id) const;
3488  global_active_cell_index() const;
3489 
3496  global_level_cell_index() const;
3497 
3511  bool
3512  direction_flag() const;
3513 
3539  unsigned int
3540  active_cell_index() const;
3541 
3549  int
3550  parent_index() const;
3551 
3558  parent() const;
3559 
3585  bool
3586  active() const;
3587 
3596  bool
3597  is_active() const;
3598 
3618  bool
3619  is_locally_owned() const;
3620 
3625  bool
3626  is_locally_owned_on_level() const;
3627 
3651  bool
3652  is_ghost() const;
3653 
3680  bool
3681  is_artificial() const;
3682 
3696  bool
3697  point_inside(const Point<spacedim> &p) const;
3698 
3707  void
3708  set_neighbor(const unsigned int i,
3709  const TriaIterator<CellAccessor<dim, spacedim>> &pointer) const;
3710 
3724  CellId
3725  id() const;
3726 
3728 
3732  double
3733  diameter(const Mapping<dim, spacedim> &mapping) const;
3734 
3743  DeclException0(ExcRefineCellNotActive);
3747  DeclException0(ExcCellFlaggedForRefinement);
3751  DeclException0(ExcCellFlaggedForCoarsening);
3752 
3753 protected:
3769  unsigned int
3770  neighbor_of_neighbor_internal(const unsigned int neighbor) const;
3771 
3777  template <int dim_, int spacedim_>
3778  bool
3779  point_inside_codim(const Point<spacedim_> &p) const;
3780 
3781 
3782 
3783 private:
3788  void
3789  set_active_cell_index(const unsigned int active_cell_index) const;
3790 
3794  void
3795  set_global_active_cell_index(const types::global_cell_index index) const;
3796 
3800  void
3801  set_global_level_cell_index(const types::global_cell_index index) const;
3802 
3806  void
3807  set_parent(const unsigned int parent_index);
3808 
3815  void
3816  set_direction_flag(const bool new_direction_flag) const;
3817 
3818  template <int, int>
3819  friend class Triangulation;
3820 
3821  template <int, int>
3823 
3824  friend struct ::internal::TriangulationImplementation::Implementation;
3825  friend struct ::internal::TriangulationImplementation::
3826  ImplementationMixedMesh;
3827 };
3828 
3829 
3830 
3831 /* ----- declaration of explicit specializations and general templates ----- */
3832 
3833 
3834 template <int structdim, int dim, int spacedim>
3835 template <typename OtherAccessor>
3837  const OtherAccessor &)
3838 {
3839  Assert(false,
3840  ExcMessage("You are attempting an illegal conversion between "
3841  "iterator/accessor types. The constructor you call "
3842  "only exists to make certain template constructs "
3843  "easier to write as dimension independent code but "
3844  "the conversion is not valid in the current context."));
3845 }
3846 
3847 
3848 
3849 template <int structdim, int dim, int spacedim>
3850 template <int structdim2, int dim2, int spacedim2>
3853 {
3854  Assert(false,
3855  ExcMessage("You are attempting an illegal conversion between "
3856  "iterator/accessor types. The constructor you call "
3857  "only exists to make certain template constructs "
3858  "easier to write as dimension independent code but "
3859  "the conversion is not valid in the current context."));
3860 }
3861 
3862 
3863 
3864 template <int dim, int spacedim>
3865 template <int structdim2, int dim2, int spacedim2>
3868 {
3869  Assert(false,
3870  ExcMessage("You are attempting an illegal conversion between "
3871  "iterator/accessor types. The constructor you call "
3872  "only exists to make certain template constructs "
3873  "easier to write as dimension independent code but "
3874  "the conversion is not valid in the current context."));
3875 }
3876 
3877 
3878 
3879 template <int structdim, int dim, int spacedim>
3880 template <int structdim2, int dim2, int spacedim2>
3883 {
3884  Assert(false,
3885  ExcMessage("You are attempting an illegal conversion between "
3886  "iterator/accessor types. The constructor you call "
3887  "only exists to make certain template constructs "
3888  "easier to write as dimension independent code but "
3889  "the conversion is not valid in the current context."));
3890 }
3891 
3892 
3893 
3894 template <int dim, int spacedim>
3895 template <int structdim2, int dim2, int spacedim2>
3898 {
3899  Assert(false,
3900  ExcMessage("You are attempting an illegal conversion between "
3901  "iterator/accessor types. The constructor you call "
3902  "only exists to make certain template constructs "
3903  "easier to write as dimension independent code but "
3904  "the conversion is not valid in the current context."));
3905 }
3906 
3907 
3908 #ifndef DOXYGEN
3909 
3910 template <>
3911 bool
3913 template <>
3914 bool
3916 template <>
3917 bool
3919 template <>
3920 bool
3922 template <>
3923 bool
3925 template <>
3926 bool
3928 // -------------------------------------------------------------------
3929 
3930 template <>
3931 void
3933 
3934 #endif // DOXYGEN
3935 
3937 
3938 // include more templates in debug and optimized mode
3939 #include "tria_accessor.templates.h"
3940 
3941 #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)
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:1466
Abstract base class for mapping classes.
Definition: mapping.h:303
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:380
unsigned int level
Definition: grid_out.cc:4370
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:70
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:379
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()