Reference documentation for deal.II version Git e694e60270 2021-02-25 22:28:05 -0500
\(\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 
33 #include <boost/container/small_vector.hpp>
35 
36 #include <utility>
37 
38 
40 
41 // Forward declarations
42 #ifndef DOXYGEN
43 template <int dim, int spacedim>
44 class Triangulation;
45 template <typename Accessor>
46 class TriaRawIterator;
47 template <typename Accessor>
48 class TriaIterator;
49 template <typename Accessor>
50 class TriaActiveIterator;
51 
52 namespace parallel
53 {
54  template <int dim, int spacedim>
55  class TriangulationBase;
56 }
57 
58 template <int dim, int spacedim>
59 class Manifold;
60 
61 template <int dim, int spacedim>
62 class Mapping;
63 #endif
64 
65 namespace internal
66 {
67  namespace TriangulationImplementation
68  {
69  class TriaObjects;
70  struct Implementation;
71  struct ImplementationMixedMesh;
72  } // namespace TriangulationImplementation
73 
74  namespace TriaAccessorImplementation
75  {
76  struct Implementation;
77 
83  template <int structdim, int dim>
85  {
86  struct type
87  {
91  type() = default;
92 
96  type(const int level)
97  {
98  Assert(level == 0, ExcInternalError());
99  (void)level; // removes -Wunused-parameter warning in optimized mode
100  }
101 
105  operator int() const
106  {
107  return 0;
108  }
109 
110  void
111  operator++() const
112  {
113  Assert(false, ExcInternalError());
114  }
115 
116  void
117  operator--() const
118  {
119  Assert(false, ExcInternalError());
120  }
121  };
122  };
123 
124 
130  template <int dim>
131  struct PresentLevelType<dim, dim>
132  {
133  using type = int;
134  };
135  } // namespace TriaAccessorImplementation
136 } // namespace internal
137 template <int structdim, int dim, int spacedim>
139 template <int dim, int spacedim>
140 class TriaAccessor<0, dim, spacedim>;
141 template <int spacedim>
142 class TriaAccessor<0, 1, spacedim>;
143 
144 // note: the file tria_accessor.templates.h is included at the end of
145 // this file. this includes a lot of templates. originally, this was
146 // only done in debug mode, but led to cyclic reduction problems and
147 // so is now on by default.
148 
149 
154 {
159  "The operation you are attempting can only be performed for "
160  "(cell, face, or edge) iterators that point to valid "
161  "objects. These objects need not necessarily be active, "
162  "i.e., have no children, but they need to be part of a "
163  "triangulation. (The objects pointed to by an iterator "
164  "may -- after coarsening -- also be objects that used "
165  "to be part of a triangulation, but are now no longer "
166  "used. Their memory location may have been retained "
167  "for re-use upon the next mesh refinement, but is "
168  "currently unused.)");
179  "The operation you are attempting can only be performed for "
180  "(cell, face, or edge) iterators that point to 'active' "
181  "objects. 'Active' objects are those that do not have "
182  "children (in the case of cells), or that are part of "
183  "an active cell (in the case of faces or edges). However, "
184  "the object on which you are trying the current "
185  "operation is not 'active' in this sense.");
192  "The operation you are attempting can only be performed for "
193  "(cell, face, or edge) iterators that have children, "
194  "but the object on which you are trying the current "
195  "operation does not have any.");
203  "The operation you are attempting can only be performed for "
204  "(cell, face, or edge) iterators that have a parent object, "
205  "but the object on which you are trying the current "
206  "operation does not have one -- i.e., it is on the "
207  "coarsest level of the triangulation.");
212  int,
213  << "You can only set the child index if the cell does not "
214  << "currently have children registered; or you can clear it. "
215  << "The given index was " << arg1
216  << " (-1 means: clear children).");
220  template <typename AccessorType>
222  AccessorType,
223  << "You tried to dereference an iterator for which this "
224  << "is not possible. More information on this iterator: "
225  << "index=" << arg1.index() << ", state="
226  << (arg1.state() == IteratorState::valid ?
227  "valid" :
228  (arg1.state() == IteratorState::past_the_end ?
229  "past_the_end" :
230  "invalid")));
235  "Iterators can only be compared if they point to the same "
236  "triangulation, or if neither of them are associated "
237  "with a triangulation.");
238  // TODO: Write documentation!
243  // TODO: Write documentation!
263  // TODO: Write documentation!
269  int,
270  << "You can only set the child index of an even numbered child."
271  << "The number of the child given was " << arg1 << ".");
272 } // namespace TriaAccessorExceptions
273 
274 
300 template <int structdim, int dim, int spacedim = dim>
302 {
303 public:
309  static const unsigned int space_dimension = spacedim;
310 
316  static const unsigned int dimension = dim;
317 
323  static const unsigned int structure_dimension = structdim;
324 
334  void
335  operator=(const TriaAccessorBase *) = delete;
336 
337 protected:
343  using AccessorData = void;
344 
348  TriaAccessorBase(const Triangulation<dim, spacedim> *parent = nullptr,
349  const int level = -1,
350  const int index = -1,
351  const AccessorData * = nullptr);
352 
357 
365  void
366  copy_from(const TriaAccessorBase &);
367 
372  operator=(const TriaAccessorBase &);
373 
384  bool
385  operator<(const TriaAccessorBase &other) const;
386 
387 protected:
391  bool
392  operator==(const TriaAccessorBase &) const;
393 
397  bool
398  operator!=(const TriaAccessorBase &) const;
399 
413  void
414  operator++();
415 
423  void
424  operator--();
433  objects() const;
434 
435 public:
441  using LocalData = void *;
442 
466  int
467  level() const;
468 
495  int
496  index() const;
497 
503  state() const;
504 
510  get_triangulation() const;
511 
515 protected:
520  typename ::internal::TriaAccessorImplementation::
521  PresentLevelType<structdim, dim>::type present_level;
522 
528 
533 
534 private:
535  template <typename Accessor>
536  friend class TriaRawIterator;
537  template <typename Accessor>
538  friend class TriaIterator;
539  template <typename Accessor>
540  friend class TriaActiveIterator;
541 };
542 
543 
544 
565 template <int structdim, int dim, int spacedim = dim>
566 class InvalidAccessor : public TriaAccessorBase<structdim, dim, spacedim>
567 {
568 public:
572  using AccessorData =
574 
582  InvalidAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
583  const int level = -1,
584  const int index = -1,
585  const AccessorData * local_data = nullptr);
586 
595 
600  template <typename OtherAccessor>
601  InvalidAccessor(const OtherAccessor &);
602 
606  void
607  copy_from(const InvalidAccessor &);
608 
612  bool
613  operator==(const InvalidAccessor &) const;
614  bool
615  operator!=(const InvalidAccessor &) const;
616 
620  void
621  operator++() const;
622  void
623  operator--() const;
624 
629  bool
630  used() const;
631 
636  bool
637  has_children() const;
638 
643  manifold_id() const;
644 
648  unsigned int
649  user_index() const;
650 
654  void
655  set_user_index(const unsigned int p) const;
656 
660  void
661  set_manifold_id(const types::manifold_id) const;
662 
667  vertex(const unsigned int i) const;
668 
673  typename ::internal::TriangulationImplementation::
674  Iterators<dim, spacedim>::line_iterator
675  line(const unsigned int i) const;
676 
681  typename ::internal::TriangulationImplementation::
682  Iterators<dim, spacedim>::quad_iterator
683  quad(const unsigned int i) const;
684 };
685 
686 
687 
705 template <int structdim, int dim, int spacedim>
706 class TriaAccessor : public TriaAccessorBase<structdim, dim, spacedim>
707 {
708 public:
712  using AccessorData =
714 
718  TriaAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
719  const int level = -1,
720  const int index = -1,
721  const AccessorData * local_data = nullptr);
722 
735  template <int structdim2, int dim2, int spacedim2>
737 
742  template <int structdim2, int dim2, int spacedim2>
744 
754  void
755  operator=(const TriaAccessor &) = delete;
756 
763  bool
764  used() const;
765 
778  vertex_iterator(const unsigned int i) const;
779 
795  unsigned int
796  vertex_index(const unsigned int i) const;
797 
836  vertex(const unsigned int i) const;
837 
841  typename ::internal::TriangulationImplementation::
842  Iterators<dim, spacedim>::line_iterator
843  line(const unsigned int i) const;
844 
851  unsigned int
852  line_index(const unsigned int i) const;
853 
857  typename ::internal::TriangulationImplementation::
858  Iterators<dim, spacedim>::quad_iterator
859  quad(const unsigned int i) const;
860 
867  unsigned int
868  quad_index(const unsigned int i) const;
891  bool
892  face_orientation(const unsigned int face) const;
893 
903  bool
904  face_flip(const unsigned int face) const;
905 
915  bool
916  face_rotation(const unsigned int face) const;
917 
928  bool
929  line_orientation(const unsigned int line) const;
944  bool
945  has_children() const;
946 
951  unsigned int
952  n_children() const;
953 
967  unsigned int
968  number_of_children() const;
969 
983  unsigned int
984  max_refinement_depth() const;
985 
990  child(const unsigned int i) const;
991 
996  unsigned int
997  child_iterator_to_index(
999 
1009  isotropic_child(const unsigned int i) const;
1010 
1015  refinement_case() const;
1016 
1022  int
1023  child_index(const unsigned int i) const;
1024 
1030  int
1031  isotropic_child_index(const unsigned int i) const;
1054  boundary_id() const;
1055 
1085  void
1086  set_boundary_id(const types::boundary_id) const;
1087 
1118  void
1119  set_all_boundary_ids(const types::boundary_id) const;
1120 
1128  bool
1129  at_boundary() const;
1130 
1140  const Manifold<dim, spacedim> &
1141  get_manifold() const;
1142 
1164  manifold_id() const;
1165 
1183  void
1184  set_manifold_id(const types::manifold_id) const;
1185 
1199  void
1200  set_all_manifold_ids(const types::manifold_id) const;
1201 
1218  bool
1219  user_flag_set() const;
1220 
1226  void
1227  set_user_flag() const;
1228 
1234  void
1235  clear_user_flag() const;
1236 
1242  void
1243  recursively_set_user_flag() const;
1244 
1250  void
1251  recursively_clear_user_flag() const;
1252 
1258  void
1259  clear_user_data() const;
1260 
1272  void
1273  set_user_pointer(void *p) const;
1274 
1280  void
1281  clear_user_pointer() const;
1282 
1298  void *
1299  user_pointer() const;
1300 
1322  void
1323  recursively_set_user_pointer(void *p) const;
1324 
1331  void
1332  recursively_clear_user_pointer() const;
1333 
1343  void
1344  set_user_index(const unsigned int p) const;
1345 
1351  void
1352  clear_user_index() const;
1353 
1365  unsigned int
1366  user_index() const;
1367 
1385  void
1386  recursively_set_user_index(const unsigned int p) const;
1387 
1396  void
1397  recursively_clear_user_index() const;
1432  double
1433  diameter() const;
1434 
1461  std::pair<Point<spacedim>, double>
1462  enclosing_ball() const;
1463 
1473  bounding_box() const;
1474 
1484  double
1485  extent_in_direction(const unsigned int axis) const;
1486 
1490  double
1491  minimum_vertex_distance() const;
1492 
1507  intermediate_point(const Point<structdim> &coordinates) const;
1508 
1532  real_to_unit_cell_affine_approximation(const Point<spacedim> &point) const;
1533 
1569  center(const bool respect_manifold = false,
1570  const bool interpolate_from_surrounding = false) const;
1571 
1590  barycenter() const;
1591 
1617  double
1618  measure() const;
1619 
1634  bool
1635  is_translation_of(
1637 
1642  reference_cell() const;
1643 
1647  unsigned int
1648  n_vertices() const;
1649 
1653  unsigned int
1654  n_lines() const;
1655 
1661  unsigned int
1662  n_faces() const;
1663 
1669  vertex_indices() const;
1670 
1676  line_indices() const;
1677 
1685  face_indices() const;
1686 
1691 private:
1696  void
1697  set_boundary_id_internal(const types::boundary_id id) const;
1698 
1706  void
1707  set_bounding_object_indices(
1708  const std::initializer_list<int> &new_indices) const;
1709 
1713  void
1714  set_bounding_object_indices(
1715  const std::initializer_list<unsigned int> &new_indices) const;
1716 
1724  void
1725  set_line_orientation(const unsigned int line, const bool orientation) const;
1726 
1737  void
1738  set_face_orientation(const unsigned int face, const bool orientation) const;
1739 
1746  void
1747  set_face_flip(const unsigned int face, const bool flip) const;
1748 
1755  void
1756  set_face_rotation(const unsigned int face, const bool rotation) const;
1757 
1761  void
1762  set_used_flag() const;
1763 
1767  void
1768  clear_used_flag() const;
1769 
1778  void
1779  set_refinement_case(const RefinementCase<structdim> &ref_case) const;
1780 
1788  void
1789  clear_refinement_case() const;
1790 
1797  void
1798  set_children(const unsigned int i, const int index) const;
1799 
1804  void
1805  clear_children() const;
1806 
1807 private:
1808  template <int, int>
1809  friend class Triangulation;
1810 
1811  friend struct ::internal::TriangulationImplementation::Implementation;
1812  friend struct ::internal::TriangulationImplementation::
1813  ImplementationMixedMesh;
1814  friend struct ::internal::TriaAccessorImplementation::Implementation;
1815 };
1816 
1817 
1818 
1837 template <int dim, int spacedim>
1838 class TriaAccessor<0, dim, spacedim>
1839 {
1840 public:
1846  static const unsigned int space_dimension = spacedim;
1847 
1853  static const unsigned int dimension = dim;
1854 
1860  static const unsigned int structure_dimension = 0;
1861 
1865  using AccessorData = void;
1866 
1872  const unsigned int vertex_index);
1873 
1879  TriaAccessor(const Triangulation<dim, spacedim> *tria = nullptr,
1880  const int level = 0,
1881  const int index = 0,
1882  const AccessorData * = nullptr);
1883 
1887  template <int structdim2, int dim2, int spacedim2>
1889 
1893  template <int structdim2, int dim2, int spacedim2>
1895 
1900  state() const;
1901 
1906  static int
1907  level();
1908 
1913  int
1914  index() const;
1915 
1921  get_triangulation() const;
1922 
1932  void
1933  operator++();
1934 
1938  void
1939  operator--();
1943  bool
1944  operator==(const TriaAccessor &) const;
1945 
1949  bool
1950  operator!=(const TriaAccessor &) const;
1951 
1979  unsigned int
1980  vertex_index(const unsigned int i = 0) const;
1981 
1987  Point<spacedim> &
1988  vertex(const unsigned int i = 0) const;
1989 
1994  typename ::internal::TriangulationImplementation::
1995  Iterators<dim, spacedim>::line_iterator static line(const unsigned int);
1996 
2000  static unsigned int
2001  line_index(const unsigned int i);
2002 
2006  static typename ::internal::TriangulationImplementation::
2007  Iterators<dim, spacedim>::quad_iterator
2008  quad(const unsigned int i);
2009 
2013  static unsigned int
2014  quad_index(const unsigned int i);
2015 
2031  double
2032  diameter() const;
2033 
2041  double
2042  extent_in_direction(const unsigned int axis) const;
2043 
2052  center(const bool respect_manifold = false,
2053  const bool interpolate_from_surrounding = false) const;
2054 
2063  double
2064  measure() const;
2079  static bool
2080  face_orientation(const unsigned int face);
2081 
2085  static bool
2086  face_flip(const unsigned int face);
2087 
2091  static bool
2092  face_rotation(const unsigned int face);
2093 
2097  static bool
2098  line_orientation(const unsigned int line);
2099 
2114  static bool
2115  has_children();
2116 
2121  static unsigned int
2122  n_children();
2123 
2128  static unsigned int
2129  number_of_children();
2130 
2134  static unsigned int
2135  max_refinement_depth();
2136 
2140  static unsigned int
2141  child_iterator_to_index(const TriaIterator<TriaAccessor<0, dim, spacedim>> &);
2142 
2147  child(const unsigned int);
2148 
2153  isotropic_child(const unsigned int);
2154 
2158  static RefinementCase<0>
2159  refinement_case();
2160 
2164  static int
2165  child_index(const unsigned int i);
2166 
2170  static int
2171  isotropic_child_index(const unsigned int i);
2179  bool
2180  used() const;
2181 
2182 protected:
2190  void
2191  copy_from(const TriaAccessor &);
2192 
2201  bool
2202  operator<(const TriaAccessor &other) const;
2203 
2208 
2212  unsigned int global_vertex_index;
2213 
2214 private:
2215  template <typename Accessor>
2216  friend class TriaRawIterator;
2217  template <typename Accessor>
2218  friend class TriaIterator;
2219  template <typename Accessor>
2220  friend class TriaActiveIterator;
2221 };
2222 
2223 
2224 
2241 template <int spacedim>
2242 class TriaAccessor<0, 1, spacedim>
2243 {
2244 public:
2250  static const unsigned int space_dimension = spacedim;
2251 
2257  static const unsigned int dimension = 1;
2258 
2264  static const unsigned int structure_dimension = 0;
2265 
2269  using AccessorData = void;
2270 
2276  {
2288  right_vertex
2289  };
2290 
2303  const VertexKind vertex_kind,
2304  const unsigned int vertex_index);
2305 
2311  TriaAccessor(const Triangulation<1, spacedim> *tria = nullptr,
2312  const int = 0,
2313  const int = 0,
2314  const AccessorData * = nullptr);
2315 
2319  template <int structdim2, int dim2, int spacedim2>
2321 
2325  template <int structdim2, int dim2, int spacedim2>
2327 
2332  void
2333  copy_from(const TriaAccessor &);
2334 
2341  state();
2342 
2347  static int
2348  level();
2349 
2354  int
2355  index() const;
2356 
2362  get_triangulation() const;
2363 
2374  void
2375  operator++() const;
2376 
2381  void
2382  operator--() const;
2386  bool
2387  operator==(const TriaAccessor &) const;
2388 
2392  bool
2393  operator!=(const TriaAccessor &) const;
2394 
2403  bool
2404  operator<(const TriaAccessor &other) const;
2405 
2432  unsigned int
2433  vertex_index(const unsigned int i = 0) const;
2434 
2440  Point<spacedim> &
2441  vertex(const unsigned int i = 0) const;
2442 
2448  center() const;
2449 
2454  typename ::internal::TriangulationImplementation::
2455  Iterators<1, spacedim>::line_iterator static line(const unsigned int);
2456 
2463  static unsigned int
2464  line_index(const unsigned int i);
2465 
2469  static typename ::internal::TriangulationImplementation::
2470  Iterators<1, spacedim>::quad_iterator
2471  quad(const unsigned int i);
2472 
2479  static unsigned int
2480  quad_index(const unsigned int i);
2481 
2491  bool
2492  at_boundary() const;
2493 
2509  boundary_id() const;
2510 
2514  const Manifold<1, spacedim> &
2515  get_manifold() const;
2516 
2524  manifold_id() const;
2525 
2526 
2537  static bool
2538  face_orientation(const unsigned int face);
2539 
2543  static bool
2544  face_flip(const unsigned int face);
2545 
2549  static bool
2550  face_rotation(const unsigned int face);
2551 
2555  static bool
2556  line_orientation(const unsigned int line);
2557 
2572  static bool
2573  has_children();
2574 
2579  static unsigned int
2580  n_children();
2581 
2586  static unsigned int
2587  number_of_children();
2588 
2592  static unsigned int
2593  max_refinement_depth();
2594 
2598  static unsigned int
2599  child_iterator_to_index(const TriaIterator<TriaAccessor<0, 1, spacedim>> &);
2600 
2605  child(const unsigned int);
2606 
2611  isotropic_child(const unsigned int);
2612 
2616  static RefinementCase<0>
2617  refinement_case();
2618 
2622  static int
2623  child_index(const unsigned int i);
2624 
2628  static int
2629  isotropic_child_index(const unsigned int i);
2662  void
2663  set_boundary_id(const types::boundary_id);
2664 
2671  void
2672  set_manifold_id(const types::manifold_id);
2673 
2685  void
2686  set_all_boundary_ids(const types::boundary_id);
2687 
2699  void
2700  set_all_manifold_ids(const types::manifold_id);
2708  bool
2709  used() const;
2710 
2715  reference_cell() const;
2716 
2720  unsigned int
2721  n_vertices() const;
2722 
2726  unsigned int
2727  n_lines() const;
2728 
2734  vertex_indices() const;
2735 
2741  line_indices() const;
2742 
2743 protected:
2748 
2754 
2758  unsigned int global_vertex_index;
2759 };
2760 
2761 
2762 
2778 template <int dim, int spacedim = dim>
2779 class CellAccessor : public TriaAccessor<dim, dim, spacedim>
2780 {
2781 public:
2786 
2791 
2802  CellAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
2803  const int level = -1,
2804  const int index = -1,
2805  const AccessorData * local_data = nullptr);
2806 
2810  CellAccessor(const TriaAccessor<dim, dim, spacedim> &cell_accessor);
2811 
2824  template <int structdim2, int dim2, int spacedim2>
2826 
2831  template <int structdim2, int dim2, int spacedim2>
2833 
2843  void
2844  operator=(const CellAccessor<dim, spacedim> &) = delete;
2845 
2862  child(const unsigned int i) const;
2863 
2867  boost::container::small_vector<TriaIterator<CellAccessor<dim, spacedim>>,
2869  child_iterators() const;
2870 
2874  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>
2875  face(const unsigned int i) const;
2876 
2881  unsigned int
2882  face_iterator_to_index(
2884 
2888  boost::container::small_vector<
2889  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>,
2891  face_iterators() const;
2892 
2902  unsigned int
2903  face_index(const unsigned int i) const;
2904 
2953  neighbor_child_on_subface(const unsigned int face_no,
2954  const unsigned int subface_no) const;
2955 
2981  neighbor(const unsigned int face_no) const;
2982 
2990  int
2991  neighbor_index(const unsigned int face_no) const;
2992 
3000  int
3001  neighbor_level(const unsigned int face_no) const;
3002 
3014  unsigned int
3015  neighbor_of_neighbor(const unsigned int face_no) const;
3016 
3027  bool
3028  neighbor_is_coarser(const unsigned int face_no) const;
3029 
3044  std::pair<unsigned int, unsigned int>
3045  neighbor_of_coarser_neighbor(const unsigned int neighbor) const;
3046 
3053  unsigned int
3054  neighbor_face_no(const unsigned int neighbor) const;
3055 
3059  static bool
3060  is_level_cell();
3061 
3075  bool
3076  has_periodic_neighbor(const unsigned int i) const;
3077 
3095  periodic_neighbor(const unsigned int i) const;
3096 
3105  neighbor_or_periodic_neighbor(const unsigned int i) const;
3106 
3122  periodic_neighbor_child_on_subface(const unsigned int face_no,
3123  const unsigned int subface_no) const;
3124 
3135  std::pair<unsigned int, unsigned int>
3136  periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const;
3137 
3143  int
3144  periodic_neighbor_index(const unsigned int i) const;
3145 
3151  int
3152  periodic_neighbor_level(const unsigned int i) const;
3153 
3168  unsigned int
3169  periodic_neighbor_of_periodic_neighbor(const unsigned int i) const;
3170 
3176  unsigned int
3177  periodic_neighbor_face_no(const unsigned int i) const;
3178 
3185  bool
3186  periodic_neighbor_is_coarser(const unsigned int i) const;
3187 
3204  bool
3205  at_boundary(const unsigned int i) const;
3206 
3215  bool
3216  at_boundary() const;
3217 
3225  bool
3226  has_boundary_lines() const;
3253  refine_flag_set() const;
3254 
3272  void
3273  set_refine_flag(const RefinementCase<dim> ref_case =
3275 
3279  void
3280  clear_refine_flag() const;
3281 
3289  bool
3290  flag_for_face_refinement(
3291  const unsigned int face_no,
3292  const RefinementCase<dim - 1> &face_refinement_case =
3294 
3300  bool
3301  flag_for_line_refinement(const unsigned int line_no) const;
3302 
3312  subface_case(const unsigned int face_no) const;
3313 
3317  bool
3318  coarsen_flag_set() const;
3319 
3324  void
3325  set_coarsen_flag() const;
3326 
3330  void
3331  clear_coarsen_flag() const;
3355  material_id() const;
3356 
3368  void
3369  set_material_id(const types::material_id new_material_id) const;
3370 
3379  void
3380  recursively_set_material_id(const types::material_id new_material_id) const;
3407  subdomain_id() const;
3408 
3424  void
3425  set_subdomain_id(const types::subdomain_id new_subdomain_id) const;
3426 
3436  level_subdomain_id() const;
3437 
3442  void
3443  set_level_subdomain_id(
3444  const types::subdomain_id new_level_subdomain_id) const;
3445 
3446 
3462  void
3463  recursively_set_subdomain_id(
3464  const types::subdomain_id new_subdomain_id) const;
3482  global_active_cell_index() const;
3483 
3490  global_level_cell_index() const;
3491 
3505  bool
3506  direction_flag() const;
3507 
3533  unsigned int
3534  active_cell_index() const;
3535 
3543  int
3544  parent_index() const;
3545 
3552  parent() const;
3553 
3579  bool
3580  active() const;
3581 
3590  bool
3591  is_active() const;
3592 
3612  bool
3613  is_locally_owned() const;
3614 
3619  bool
3620  is_locally_owned_on_level() const;
3621 
3645  bool
3646  is_ghost() const;
3647 
3674  bool
3675  is_artificial() const;
3676 
3690  bool
3691  point_inside(const Point<spacedim> &p) const;
3692 
3701  void
3702  set_neighbor(const unsigned int i,
3703  const TriaIterator<CellAccessor<dim, spacedim>> &pointer) const;
3704 
3718  CellId
3719  id() const;
3720 
3722 
3726  double
3727  diameter(const Mapping<dim, spacedim> &mapping) const;
3728 
3737  DeclException0(ExcRefineCellNotActive);
3741  DeclException0(ExcCellFlaggedForRefinement);
3745  DeclException0(ExcCellFlaggedForCoarsening);
3746 
3747 protected:
3763  unsigned int
3764  neighbor_of_neighbor_internal(const unsigned int neighbor) const;
3765 
3771  template <int dim_, int spacedim_>
3772  bool
3773  point_inside_codim(const Point<spacedim_> &p) const;
3774 
3775 
3776 
3777 private:
3782  void
3783  set_active_cell_index(const unsigned int active_cell_index) const;
3784 
3788  void
3789  set_global_active_cell_index(const types::global_cell_index index) const;
3790 
3794  void
3795  set_global_level_cell_index(const types::global_cell_index index) const;
3796 
3800  void
3801  set_parent(const unsigned int parent_index);
3802 
3809  void
3810  set_direction_flag(const bool new_direction_flag) const;
3811 
3812  template <int, int>
3813  friend class Triangulation;
3814 
3815  template <int, int>
3817 
3818  friend struct ::internal::TriangulationImplementation::Implementation;
3819  friend struct ::internal::TriangulationImplementation::
3820  ImplementationMixedMesh;
3821 };
3822 
3823 
3824 
3825 /* ----- declaration of explicit specializations and general templates ----- */
3826 
3827 
3828 template <int structdim, int dim, int spacedim>
3829 template <typename OtherAccessor>
3831  const OtherAccessor &)
3832 {
3833  Assert(false,
3834  ExcMessage("You are attempting an illegal conversion between "
3835  "iterator/accessor types. The constructor you call "
3836  "only exists to make certain template constructs "
3837  "easier to write as dimension independent code but "
3838  "the conversion is not valid in the current context."));
3839 }
3840 
3841 
3842 
3843 template <int structdim, int dim, int spacedim>
3844 template <int structdim2, int dim2, int spacedim2>
3847 {
3848  Assert(false,
3849  ExcMessage("You are attempting an illegal conversion between "
3850  "iterator/accessor types. The constructor you call "
3851  "only exists to make certain template constructs "
3852  "easier to write as dimension independent code but "
3853  "the conversion is not valid in the current context."));
3854 }
3855 
3856 
3857 
3858 template <int dim, int spacedim>
3859 template <int structdim2, int dim2, int spacedim2>
3862 {
3863  Assert(false,
3864  ExcMessage("You are attempting an illegal conversion between "
3865  "iterator/accessor types. The constructor you call "
3866  "only exists to make certain template constructs "
3867  "easier to write as dimension independent code but "
3868  "the conversion is not valid in the current context."));
3869 }
3870 
3871 
3872 
3873 template <int structdim, int dim, int spacedim>
3874 template <int structdim2, int dim2, int spacedim2>
3877 {
3878  Assert(false,
3879  ExcMessage("You are attempting an illegal conversion between "
3880  "iterator/accessor types. The constructor you call "
3881  "only exists to make certain template constructs "
3882  "easier to write as dimension independent code but "
3883  "the conversion is not valid in the current context."));
3884 }
3885 
3886 
3887 
3888 template <int dim, int spacedim>
3889 template <int structdim2, int dim2, int spacedim2>
3892 {
3893  Assert(false,
3894  ExcMessage("You are attempting an illegal conversion between "
3895  "iterator/accessor types. The constructor you call "
3896  "only exists to make certain template constructs "
3897  "easier to write as dimension independent code but "
3898  "the conversion is not valid in the current context."));
3899 }
3900 
3901 
3902 #ifndef DOXYGEN
3903 
3904 template <>
3905 bool
3907 template <>
3908 bool
3910 template <>
3911 bool
3913 template <>
3914 bool
3916 template <>
3917 bool
3919 template <>
3920 bool
3922 // -------------------------------------------------------------------
3923 
3924 template <>
3925 void
3927 
3928 #endif // DOXYGEN
3929 
3931 
3932 // include more templates in debug and optimized mode
3933 #include "tria_accessor.templates.h"
3934 
3935 #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:80
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()
#define DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
Definition: config.h:407
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
void reference_cell(const ReferenceCell &reference_cell, Triangulation< dim, spacedim > &tria)
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:394
unsigned int level
Definition: grid_out.cc:4578
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_ENABLE_EXTRA_DIAGNOSTICS
Definition: config.h:444
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:393
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:156
unsigned int boundary_id
Definition: types.h:129
typename ::internal::TriaAccessorImplementation::PresentLevelType< structdim, dim >::type present_level
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcCantCompareIterators()