Reference documentation for deal.II version Git 87130829be 2021-04-14 12:03:57 -0400
\(\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 
930  bool
931  line_orientation(const unsigned int line) const;
946  bool
947  has_children() const;
948 
953  unsigned int
954  n_children() const;
955 
969  unsigned int
970  number_of_children() const;
971 
985  unsigned int
986  max_refinement_depth() const;
987 
992  child(const unsigned int i) const;
993 
998  unsigned int
999  child_iterator_to_index(
1001 
1011  isotropic_child(const unsigned int i) const;
1012 
1017  refinement_case() const;
1018 
1024  int
1025  child_index(const unsigned int i) const;
1026 
1032  int
1033  isotropic_child_index(const unsigned int i) const;
1056  boundary_id() const;
1057 
1087  void
1088  set_boundary_id(const types::boundary_id) const;
1089 
1120  void
1121  set_all_boundary_ids(const types::boundary_id) const;
1122 
1130  bool
1131  at_boundary() const;
1132 
1142  const Manifold<dim, spacedim> &
1143  get_manifold() const;
1144 
1166  manifold_id() const;
1167 
1185  void
1186  set_manifold_id(const types::manifold_id) const;
1187 
1201  void
1202  set_all_manifold_ids(const types::manifold_id) const;
1203 
1220  bool
1221  user_flag_set() const;
1222 
1228  void
1229  set_user_flag() const;
1230 
1236  void
1237  clear_user_flag() const;
1238 
1244  void
1245  recursively_set_user_flag() const;
1246 
1252  void
1253  recursively_clear_user_flag() const;
1254 
1260  void
1261  clear_user_data() const;
1262 
1274  void
1275  set_user_pointer(void *p) const;
1276 
1282  void
1283  clear_user_pointer() const;
1284 
1300  void *
1301  user_pointer() const;
1302 
1324  void
1325  recursively_set_user_pointer(void *p) const;
1326 
1333  void
1334  recursively_clear_user_pointer() const;
1335 
1345  void
1346  set_user_index(const unsigned int p) const;
1347 
1353  void
1354  clear_user_index() const;
1355 
1367  unsigned int
1368  user_index() const;
1369 
1387  void
1388  recursively_set_user_index(const unsigned int p) const;
1389 
1398  void
1399  recursively_clear_user_index() const;
1434  double
1435  diameter() const;
1436 
1463  std::pair<Point<spacedim>, double>
1464  enclosing_ball() const;
1465 
1475  bounding_box() const;
1476 
1486  double
1487  extent_in_direction(const unsigned int axis) const;
1488 
1492  double
1493  minimum_vertex_distance() const;
1494 
1509  intermediate_point(const Point<structdim> &coordinates) const;
1510 
1534  real_to_unit_cell_affine_approximation(const Point<spacedim> &point) const;
1535 
1571  center(const bool respect_manifold = false,
1572  const bool interpolate_from_surrounding = false) const;
1573 
1592  barycenter() const;
1593 
1619  double
1620  measure() const;
1621 
1636  bool
1637  is_translation_of(
1639 
1644  reference_cell() const;
1645 
1649  unsigned int
1650  n_vertices() const;
1651 
1655  unsigned int
1656  n_lines() const;
1657 
1663  unsigned int
1664  n_faces() const;
1665 
1671  vertex_indices() const;
1672 
1678  line_indices() const;
1679 
1687  face_indices() const;
1688 
1693 private:
1698  void
1699  set_boundary_id_internal(const types::boundary_id id) const;
1700 
1708  void
1709  set_bounding_object_indices(
1710  const std::initializer_list<int> &new_indices) const;
1711 
1715  void
1716  set_bounding_object_indices(
1717  const std::initializer_list<unsigned int> &new_indices) const;
1718 
1726  void
1727  set_line_orientation(const unsigned int line, const bool orientation) const;
1728 
1739  void
1740  set_face_orientation(const unsigned int face, const bool orientation) const;
1741 
1748  void
1749  set_face_flip(const unsigned int face, const bool flip) const;
1750 
1757  void
1758  set_face_rotation(const unsigned int face, const bool rotation) const;
1759 
1763  void
1764  set_used_flag() const;
1765 
1769  void
1770  clear_used_flag() const;
1771 
1780  void
1781  set_refinement_case(const RefinementCase<structdim> &ref_case) const;
1782 
1790  void
1791  clear_refinement_case() const;
1792 
1799  void
1800  set_children(const unsigned int i, const int index) const;
1801 
1806  void
1807  clear_children() const;
1808 
1809 private:
1810  template <int, int>
1811  friend class Triangulation;
1812 
1813  friend struct ::internal::TriangulationImplementation::Implementation;
1814  friend struct ::internal::TriangulationImplementation::
1815  ImplementationMixedMesh;
1816  friend struct ::internal::TriaAccessorImplementation::Implementation;
1817 };
1818 
1819 
1820 
1839 template <int dim, int spacedim>
1840 class TriaAccessor<0, dim, spacedim>
1841 {
1842 public:
1848  static const unsigned int space_dimension = spacedim;
1849 
1855  static const unsigned int dimension = dim;
1856 
1862  static const unsigned int structure_dimension = 0;
1863 
1867  using AccessorData = void;
1868 
1874  const unsigned int vertex_index);
1875 
1881  TriaAccessor(const Triangulation<dim, spacedim> *tria = nullptr,
1882  const int level = 0,
1883  const int index = 0,
1884  const AccessorData * = nullptr);
1885 
1889  template <int structdim2, int dim2, int spacedim2>
1891 
1895  template <int structdim2, int dim2, int spacedim2>
1897 
1902  state() const;
1903 
1908  static int
1909  level();
1910 
1915  int
1916  index() const;
1917 
1923  get_triangulation() const;
1924 
1934  void
1935  operator++();
1936 
1940  void
1941  operator--();
1945  bool
1946  operator==(const TriaAccessor &) const;
1947 
1951  bool
1952  operator!=(const TriaAccessor &) const;
1953 
1981  unsigned int
1982  vertex_index(const unsigned int i = 0) const;
1983 
1989  Point<spacedim> &
1990  vertex(const unsigned int i = 0) const;
1991 
1996  typename ::internal::TriangulationImplementation::
1997  Iterators<dim, spacedim>::line_iterator static line(const unsigned int);
1998 
2002  static unsigned int
2003  line_index(const unsigned int i);
2004 
2008  static typename ::internal::TriangulationImplementation::
2009  Iterators<dim, spacedim>::quad_iterator
2010  quad(const unsigned int i);
2011 
2015  static unsigned int
2016  quad_index(const unsigned int i);
2017 
2033  double
2034  diameter() const;
2035 
2043  double
2044  extent_in_direction(const unsigned int axis) const;
2045 
2054  center(const bool respect_manifold = false,
2055  const bool interpolate_from_surrounding = false) const;
2056 
2065  double
2066  measure() const;
2081  static bool
2082  face_orientation(const unsigned int face);
2083 
2087  static bool
2088  face_flip(const unsigned int face);
2089 
2093  static bool
2094  face_rotation(const unsigned int face);
2095 
2099  static bool
2100  line_orientation(const unsigned int line);
2101 
2116  static bool
2117  has_children();
2118 
2123  static unsigned int
2124  n_children();
2125 
2130  static unsigned int
2131  number_of_children();
2132 
2136  static unsigned int
2137  max_refinement_depth();
2138 
2142  static unsigned int
2143  child_iterator_to_index(const TriaIterator<TriaAccessor<0, dim, spacedim>> &);
2144 
2149  child(const unsigned int);
2150 
2155  isotropic_child(const unsigned int);
2156 
2160  static RefinementCase<0>
2161  refinement_case();
2162 
2166  static int
2167  child_index(const unsigned int i);
2168 
2172  static int
2173  isotropic_child_index(const unsigned int i);
2181  bool
2182  used() const;
2183 
2184 protected:
2192  void
2193  copy_from(const TriaAccessor &);
2194 
2203  bool
2204  operator<(const TriaAccessor &other) const;
2205 
2210 
2214  unsigned int global_vertex_index;
2215 
2216 private:
2217  template <typename Accessor>
2218  friend class TriaRawIterator;
2219  template <typename Accessor>
2220  friend class TriaIterator;
2221  template <typename Accessor>
2222  friend class TriaActiveIterator;
2223 };
2224 
2225 
2226 
2243 template <int spacedim>
2244 class TriaAccessor<0, 1, spacedim>
2245 {
2246 public:
2252  static const unsigned int space_dimension = spacedim;
2253 
2259  static const unsigned int dimension = 1;
2260 
2266  static const unsigned int structure_dimension = 0;
2267 
2271  using AccessorData = void;
2272 
2278  {
2290  right_vertex
2291  };
2292 
2305  const VertexKind vertex_kind,
2306  const unsigned int vertex_index);
2307 
2313  TriaAccessor(const Triangulation<1, spacedim> *tria = nullptr,
2314  const int = 0,
2315  const int = 0,
2316  const AccessorData * = nullptr);
2317 
2321  template <int structdim2, int dim2, int spacedim2>
2323 
2327  template <int structdim2, int dim2, int spacedim2>
2329 
2334  void
2335  copy_from(const TriaAccessor &);
2336 
2343  state();
2344 
2349  static int
2350  level();
2351 
2356  int
2357  index() const;
2358 
2364  get_triangulation() const;
2365 
2376  void
2377  operator++() const;
2378 
2383  void
2384  operator--() const;
2388  bool
2389  operator==(const TriaAccessor &) const;
2390 
2394  bool
2395  operator!=(const TriaAccessor &) const;
2396 
2405  bool
2406  operator<(const TriaAccessor &other) const;
2407 
2434  unsigned int
2435  vertex_index(const unsigned int i = 0) const;
2436 
2442  Point<spacedim> &
2443  vertex(const unsigned int i = 0) const;
2444 
2450  center() const;
2451 
2456  typename ::internal::TriangulationImplementation::
2457  Iterators<1, spacedim>::line_iterator static line(const unsigned int);
2458 
2465  static unsigned int
2466  line_index(const unsigned int i);
2467 
2471  static typename ::internal::TriangulationImplementation::
2472  Iterators<1, spacedim>::quad_iterator
2473  quad(const unsigned int i);
2474 
2481  static unsigned int
2482  quad_index(const unsigned int i);
2483 
2493  bool
2494  at_boundary() const;
2495 
2511  boundary_id() const;
2512 
2516  const Manifold<1, spacedim> &
2517  get_manifold() const;
2518 
2526  manifold_id() const;
2527 
2528 
2539  static bool
2540  face_orientation(const unsigned int face);
2541 
2545  static bool
2546  face_flip(const unsigned int face);
2547 
2551  static bool
2552  face_rotation(const unsigned int face);
2553 
2557  static bool
2558  line_orientation(const unsigned int line);
2559 
2574  static bool
2575  has_children();
2576 
2581  static unsigned int
2582  n_children();
2583 
2588  static unsigned int
2589  number_of_children();
2590 
2594  static unsigned int
2595  max_refinement_depth();
2596 
2600  static unsigned int
2601  child_iterator_to_index(const TriaIterator<TriaAccessor<0, 1, spacedim>> &);
2602 
2607  child(const unsigned int);
2608 
2613  isotropic_child(const unsigned int);
2614 
2618  static RefinementCase<0>
2619  refinement_case();
2620 
2624  static int
2625  child_index(const unsigned int i);
2626 
2630  static int
2631  isotropic_child_index(const unsigned int i);
2664  void
2665  set_boundary_id(const types::boundary_id);
2666 
2673  void
2674  set_manifold_id(const types::manifold_id);
2675 
2687  void
2688  set_all_boundary_ids(const types::boundary_id);
2689 
2701  void
2702  set_all_manifold_ids(const types::manifold_id);
2710  bool
2711  used() const;
2712 
2717  reference_cell() const;
2718 
2722  unsigned int
2723  n_vertices() const;
2724 
2728  unsigned int
2729  n_lines() const;
2730 
2736  vertex_indices() const;
2737 
2743  line_indices() const;
2744 
2745 protected:
2750 
2756 
2760  unsigned int global_vertex_index;
2761 };
2762 
2763 
2764 
2780 template <int dim, int spacedim = dim>
2781 class CellAccessor : public TriaAccessor<dim, dim, spacedim>
2782 {
2783 public:
2788 
2793 
2804  CellAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
2805  const int level = -1,
2806  const int index = -1,
2807  const AccessorData * local_data = nullptr);
2808 
2812  CellAccessor(const TriaAccessor<dim, dim, spacedim> &cell_accessor);
2813 
2826  template <int structdim2, int dim2, int spacedim2>
2828 
2833  template <int structdim2, int dim2, int spacedim2>
2835 
2845  void
2846  operator=(const CellAccessor<dim, spacedim> &) = delete;
2847 
2864  child(const unsigned int i) const;
2865 
2869  boost::container::small_vector<TriaIterator<CellAccessor<dim, spacedim>>,
2871  child_iterators() const;
2872 
2876  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>
2877  face(const unsigned int i) const;
2878 
2883  unsigned int
2884  face_iterator_to_index(
2886 
2890  boost::container::small_vector<
2891  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>,
2893  face_iterators() const;
2894 
2904  unsigned int
2905  face_index(const unsigned int i) const;
2906 
2955  neighbor_child_on_subface(const unsigned int face_no,
2956  const unsigned int subface_no) const;
2957 
2983  neighbor(const unsigned int face_no) const;
2984 
2992  int
2993  neighbor_index(const unsigned int face_no) const;
2994 
3002  int
3003  neighbor_level(const unsigned int face_no) const;
3004 
3016  unsigned int
3017  neighbor_of_neighbor(const unsigned int face_no) const;
3018 
3029  bool
3030  neighbor_is_coarser(const unsigned int face_no) const;
3031 
3046  std::pair<unsigned int, unsigned int>
3047  neighbor_of_coarser_neighbor(const unsigned int neighbor) const;
3048 
3055  unsigned int
3056  neighbor_face_no(const unsigned int neighbor) const;
3057 
3061  static bool
3062  is_level_cell();
3063 
3077  bool
3078  has_periodic_neighbor(const unsigned int i) const;
3079 
3097  periodic_neighbor(const unsigned int i) const;
3098 
3107  neighbor_or_periodic_neighbor(const unsigned int i) const;
3108 
3124  periodic_neighbor_child_on_subface(const unsigned int face_no,
3125  const unsigned int subface_no) const;
3126 
3137  std::pair<unsigned int, unsigned int>
3138  periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const;
3139 
3145  int
3146  periodic_neighbor_index(const unsigned int i) const;
3147 
3153  int
3154  periodic_neighbor_level(const unsigned int i) const;
3155 
3170  unsigned int
3171  periodic_neighbor_of_periodic_neighbor(const unsigned int i) const;
3172 
3178  unsigned int
3179  periodic_neighbor_face_no(const unsigned int i) const;
3180 
3187  bool
3188  periodic_neighbor_is_coarser(const unsigned int i) const;
3189 
3206  bool
3207  at_boundary(const unsigned int i) const;
3208 
3217  bool
3218  at_boundary() const;
3219 
3227  bool
3228  has_boundary_lines() const;
3255  refine_flag_set() const;
3256 
3274  void
3275  set_refine_flag(const RefinementCase<dim> ref_case =
3277 
3281  void
3282  clear_refine_flag() const;
3283 
3291  bool
3292  flag_for_face_refinement(
3293  const unsigned int face_no,
3294  const RefinementCase<dim - 1> &face_refinement_case =
3296 
3302  bool
3303  flag_for_line_refinement(const unsigned int line_no) const;
3304 
3314  subface_case(const unsigned int face_no) const;
3315 
3319  bool
3320  coarsen_flag_set() const;
3321 
3326  void
3327  set_coarsen_flag() const;
3328 
3332  void
3333  clear_coarsen_flag() const;
3357  material_id() const;
3358 
3370  void
3371  set_material_id(const types::material_id new_material_id) const;
3372 
3381  void
3382  recursively_set_material_id(const types::material_id new_material_id) const;
3409  subdomain_id() const;
3410 
3426  void
3427  set_subdomain_id(const types::subdomain_id new_subdomain_id) const;
3428 
3438  level_subdomain_id() const;
3439 
3444  void
3445  set_level_subdomain_id(
3446  const types::subdomain_id new_level_subdomain_id) const;
3447 
3448 
3464  void
3465  recursively_set_subdomain_id(
3466  const types::subdomain_id new_subdomain_id) const;
3484  global_active_cell_index() const;
3485 
3492  global_level_cell_index() const;
3493 
3507  bool
3508  direction_flag() const;
3509 
3535  unsigned int
3536  active_cell_index() const;
3537 
3545  int
3546  parent_index() const;
3547 
3554  parent() const;
3555 
3581  bool
3582  active() const;
3583 
3592  bool
3593  is_active() const;
3594 
3614  bool
3615  is_locally_owned() const;
3616 
3621  bool
3622  is_locally_owned_on_level() const;
3623 
3647  bool
3648  is_ghost() const;
3649 
3676  bool
3677  is_artificial() const;
3678 
3692  bool
3693  point_inside(const Point<spacedim> &p) const;
3694 
3703  void
3704  set_neighbor(const unsigned int i,
3705  const TriaIterator<CellAccessor<dim, spacedim>> &pointer) const;
3706 
3720  CellId
3721  id() const;
3722 
3724 
3728  double
3729  diameter(const Mapping<dim, spacedim> &mapping) const;
3730 
3739  DeclException0(ExcRefineCellNotActive);
3743  DeclException0(ExcCellFlaggedForRefinement);
3747  DeclException0(ExcCellFlaggedForCoarsening);
3748 
3749 protected:
3765  unsigned int
3766  neighbor_of_neighbor_internal(const unsigned int neighbor) const;
3767 
3773  template <int dim_, int spacedim_>
3774  bool
3775  point_inside_codim(const Point<spacedim_> &p) const;
3776 
3777 
3778 
3779 private:
3784  void
3785  set_active_cell_index(const unsigned int active_cell_index) const;
3786 
3790  void
3791  set_global_active_cell_index(const types::global_cell_index index) const;
3792 
3796  void
3797  set_global_level_cell_index(const types::global_cell_index index) const;
3798 
3802  void
3803  set_parent(const unsigned int parent_index);
3804 
3811  void
3812  set_direction_flag(const bool new_direction_flag) const;
3813 
3814  template <int, int>
3815  friend class Triangulation;
3816 
3817  template <int, int>
3819 
3820  friend struct ::internal::TriangulationImplementation::Implementation;
3821  friend struct ::internal::TriangulationImplementation::
3822  ImplementationMixedMesh;
3823 };
3824 
3825 
3826 
3827 /* ----- declaration of explicit specializations and general templates ----- */
3828 
3829 
3830 template <int structdim, int dim, int spacedim>
3831 template <typename OtherAccessor>
3833  const OtherAccessor &)
3834 {
3835  Assert(false,
3836  ExcMessage("You are attempting an illegal conversion between "
3837  "iterator/accessor types. The constructor you call "
3838  "only exists to make certain template constructs "
3839  "easier to write as dimension independent code but "
3840  "the conversion is not valid in the current context."));
3841 }
3842 
3843 
3844 
3845 template <int structdim, int dim, int spacedim>
3846 template <int structdim2, int dim2, int spacedim2>
3849 {
3850  Assert(false,
3851  ExcMessage("You are attempting an illegal conversion between "
3852  "iterator/accessor types. The constructor you call "
3853  "only exists to make certain template constructs "
3854  "easier to write as dimension independent code but "
3855  "the conversion is not valid in the current context."));
3856 }
3857 
3858 
3859 
3860 template <int dim, int spacedim>
3861 template <int structdim2, int dim2, int spacedim2>
3864 {
3865  Assert(false,
3866  ExcMessage("You are attempting an illegal conversion between "
3867  "iterator/accessor types. The constructor you call "
3868  "only exists to make certain template constructs "
3869  "easier to write as dimension independent code but "
3870  "the conversion is not valid in the current context."));
3871 }
3872 
3873 
3874 
3875 template <int structdim, int dim, int spacedim>
3876 template <int structdim2, int dim2, int spacedim2>
3879 {
3880  Assert(false,
3881  ExcMessage("You are attempting an illegal conversion between "
3882  "iterator/accessor types. The constructor you call "
3883  "only exists to make certain template constructs "
3884  "easier to write as dimension independent code but "
3885  "the conversion is not valid in the current context."));
3886 }
3887 
3888 
3889 
3890 template <int dim, int spacedim>
3891 template <int structdim2, int dim2, int spacedim2>
3894 {
3895  Assert(false,
3896  ExcMessage("You are attempting an illegal conversion between "
3897  "iterator/accessor types. The constructor you call "
3898  "only exists to make certain template constructs "
3899  "easier to write as dimension independent code but "
3900  "the conversion is not valid in the current context."));
3901 }
3902 
3903 
3904 #ifndef DOXYGEN
3905 
3906 template <>
3907 bool
3909 template <>
3910 bool
3912 template <>
3913 bool
3915 template <>
3916 bool
3918 template <>
3919 bool
3921 template <>
3922 bool
3924 // -------------------------------------------------------------------
3925 
3926 template <>
3927 void
3929 
3930 #endif // DOXYGEN
3931 
3933 
3934 // include more templates in debug and optimized mode
3935 #include "tria_accessor.templates.h"
3936 
3937 #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:81
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:409
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:396
unsigned int level
Definition: grid_out.cc:4590
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:446
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:395
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:158
unsigned int boundary_id
Definition: types.h:129
typename ::internal::TriaAccessorImplementation::PresentLevelType< structdim, dim >::type present_level
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcCantCompareIterators()