Reference documentation for deal.II version Git d3aed38b93 2021-10-28 13:33:27 +0200
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
tria_accessor.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2021 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 
727  TriaAccessor(const TriaAccessor &) = default;
728 
732  TriaAccessor(TriaAccessor &&) = default; // NOLINT
733 
746  template <int structdim2, int dim2, int spacedim2>
748 
753  template <int structdim2, int dim2, int spacedim2>
755 
765  TriaAccessor &
766  operator=(const TriaAccessor &) = delete;
767 
771  TriaAccessor &
772  operator=(TriaAccessor &&) = default; // NOLINT
773 
777  ~TriaAccessor() = default;
778 
785  bool
786  used() const;
787 
800  vertex_iterator(const unsigned int i) const;
801 
817  unsigned int
818  vertex_index(const unsigned int i) const;
819 
858  vertex(const unsigned int i) const;
859 
863  typename ::internal::TriangulationImplementation::
864  Iterators<dim, spacedim>::line_iterator
865  line(const unsigned int i) const;
866 
873  unsigned int
874  line_index(const unsigned int i) const;
875 
879  typename ::internal::TriangulationImplementation::
880  Iterators<dim, spacedim>::quad_iterator
881  quad(const unsigned int i) const;
882 
889  unsigned int
890  quad_index(const unsigned int i) const;
906  unsigned char
907  combined_face_orientation(const unsigned int face) const;
908 
920  bool
921  face_orientation(const unsigned int face) const;
922 
932  bool
933  face_flip(const unsigned int face) const;
934 
944  bool
945  face_rotation(const unsigned int face) const;
946 
959  bool
960  line_orientation(const unsigned int line) const;
975  bool
976  has_children() const;
977 
982  unsigned int
983  n_children() const;
984 
989  unsigned int
990  number_of_children() const;
991 
1005  unsigned int
1006  n_active_descendants() const;
1007 
1021  unsigned int
1022  max_refinement_depth() const;
1023 
1028  child(const unsigned int i) const;
1029 
1034  unsigned int
1035  child_iterator_to_index(
1037 
1047  isotropic_child(const unsigned int i) const;
1048 
1053  refinement_case() const;
1054 
1060  int
1061  child_index(const unsigned int i) const;
1062 
1068  int
1069  isotropic_child_index(const unsigned int i) const;
1092  boundary_id() const;
1093 
1123  void
1124  set_boundary_id(const types::boundary_id) const;
1125 
1156  void
1157  set_all_boundary_ids(const types::boundary_id) const;
1158 
1166  bool
1167  at_boundary() const;
1168 
1178  const Manifold<dim, spacedim> &
1179  get_manifold() const;
1180 
1202  manifold_id() const;
1203 
1221  void
1222  set_manifold_id(const types::manifold_id) const;
1223 
1237  void
1238  set_all_manifold_ids(const types::manifold_id) const;
1239 
1256  bool
1257  user_flag_set() const;
1258 
1264  void
1265  set_user_flag() const;
1266 
1272  void
1273  clear_user_flag() const;
1274 
1280  void
1281  recursively_set_user_flag() const;
1282 
1288  void
1289  recursively_clear_user_flag() const;
1290 
1296  void
1297  clear_user_data() const;
1298 
1310  void
1311  set_user_pointer(void *p) const;
1312 
1318  void
1319  clear_user_pointer() const;
1320 
1336  void *
1337  user_pointer() const;
1338 
1360  void
1361  recursively_set_user_pointer(void *p) const;
1362 
1369  void
1370  recursively_clear_user_pointer() const;
1371 
1381  void
1382  set_user_index(const unsigned int p) const;
1383 
1389  void
1390  clear_user_index() const;
1391 
1403  unsigned int
1404  user_index() const;
1405 
1423  void
1424  recursively_set_user_index(const unsigned int p) const;
1425 
1434  void
1435  recursively_clear_user_index() const;
1470  double
1471  diameter() const;
1472 
1499  std::pair<Point<spacedim>, double>
1500  enclosing_ball() const;
1501 
1511  bounding_box() const;
1512 
1522  double
1523  extent_in_direction(const unsigned int axis) const;
1524 
1528  double
1529  minimum_vertex_distance() const;
1530 
1545  intermediate_point(const Point<structdim> &coordinates) const;
1546 
1570  real_to_unit_cell_affine_approximation(const Point<spacedim> &point) const;
1571 
1607  center(const bool respect_manifold = false,
1608  const bool interpolate_from_surrounding = false) const;
1609 
1628  barycenter() const;
1629 
1655  double
1656  measure() const;
1657 
1672  bool
1673  is_translation_of(
1675 
1680  reference_cell() const;
1681 
1685  unsigned int
1686  n_vertices() const;
1687 
1691  unsigned int
1692  n_lines() const;
1693 
1699  unsigned int
1700  n_faces() const;
1701 
1707  vertex_indices() const;
1708 
1714  line_indices() const;
1715 
1723  face_indices() const;
1724 
1729 private:
1734  void
1735  set_boundary_id_internal(const types::boundary_id id) const;
1736 
1744  void
1745  set_bounding_object_indices(
1746  const std::initializer_list<int> &new_indices) const;
1747 
1751  void
1752  set_bounding_object_indices(
1753  const std::initializer_list<unsigned int> &new_indices) const;
1754 
1762  void
1763  set_line_orientation(const unsigned int line, const bool orientation) const;
1764 
1775  void
1776  set_face_orientation(const unsigned int face, const bool orientation) const;
1777 
1784  void
1785  set_face_flip(const unsigned int face, const bool flip) const;
1786 
1793  void
1794  set_face_rotation(const unsigned int face, const bool rotation) const;
1795 
1799  void
1800  set_used_flag() const;
1801 
1805  void
1806  clear_used_flag() const;
1807 
1816  void
1817  set_refinement_case(const RefinementCase<structdim> &ref_case) const;
1818 
1826  void
1827  clear_refinement_case() const;
1828 
1835  void
1836  set_children(const unsigned int i, const int index) const;
1837 
1842  void
1843  clear_children() const;
1844 
1845 private:
1846  template <int, int>
1847  friend class Triangulation;
1848 
1849  friend struct ::internal::TriangulationImplementation::Implementation;
1850  friend struct ::internal::TriangulationImplementation::
1851  ImplementationMixedMesh;
1852  friend struct ::internal::TriaAccessorImplementation::Implementation;
1853 };
1854 
1855 
1856 
1875 template <int dim, int spacedim>
1876 class TriaAccessor<0, dim, spacedim>
1877 {
1878 public:
1884  static const unsigned int space_dimension = spacedim;
1885 
1891  static const unsigned int dimension = dim;
1892 
1898  static const unsigned int structure_dimension = 0;
1899 
1903  using AccessorData = void;
1904 
1910  const unsigned int vertex_index);
1911 
1918  const int level = 0,
1919  const int index = 0,
1920  const AccessorData * = nullptr);
1921 
1925  template <int structdim2, int dim2, int spacedim2>
1927 
1931  template <int structdim2, int dim2, int spacedim2>
1933 
1938  state() const;
1939 
1944  static int
1945  level();
1946 
1951  int
1952  index() const;
1953 
1959  get_triangulation() const;
1960 
1970  void
1971  operator++();
1972 
1976  void
1977  operator--();
1981  bool
1982  operator==(const TriaAccessor &) const;
1983 
1987  bool
1988  operator!=(const TriaAccessor &) const;
1989 
2017  unsigned int
2018  vertex_index(const unsigned int i = 0) const;
2019 
2025  Point<spacedim> &
2026  vertex(const unsigned int i = 0) const;
2027 
2032  typename ::internal::TriangulationImplementation::
2033  Iterators<dim, spacedim>::line_iterator static line(const unsigned int);
2034 
2038  static unsigned int
2039  line_index(const unsigned int i);
2040 
2044  static typename ::internal::TriangulationImplementation::
2045  Iterators<dim, spacedim>::quad_iterator
2046  quad(const unsigned int i);
2047 
2051  static unsigned int
2052  quad_index(const unsigned int i);
2053 
2069  double
2070  diameter() const;
2071 
2079  double
2080  extent_in_direction(const unsigned int axis) const;
2081 
2090  center(const bool respect_manifold = false,
2091  const bool interpolate_from_surrounding = false) const;
2092 
2101  double
2102  measure() const;
2117  static unsigned char
2118  combined_face_orientation(const unsigned int face);
2119 
2123  static bool
2124  face_orientation(const unsigned int face);
2125 
2129  static bool
2130  face_flip(const unsigned int face);
2131 
2135  static bool
2136  face_rotation(const unsigned int face);
2137 
2141  static bool
2142  line_orientation(const unsigned int line);
2143 
2158  static bool
2159  has_children();
2160 
2165  static unsigned int
2166  n_children();
2167 
2172  static unsigned int
2173  n_active_descendants();
2174 
2179  static unsigned int
2180  number_of_children();
2181 
2182 
2186  static unsigned int
2187  max_refinement_depth();
2188 
2192  static unsigned int
2193  child_iterator_to_index(const TriaIterator<TriaAccessor<0, dim, spacedim>> &);
2194 
2199  child(const unsigned int);
2200 
2205  isotropic_child(const unsigned int);
2206 
2210  static RefinementCase<0>
2211  refinement_case();
2212 
2216  static int
2217  child_index(const unsigned int i);
2218 
2222  static int
2223  isotropic_child_index(const unsigned int i);
2231  bool
2232  used() const;
2233 
2234 protected:
2242  void
2243  copy_from(const TriaAccessor &);
2244 
2253  bool
2254  operator<(const TriaAccessor &other) const;
2255 
2260 
2264  unsigned int global_vertex_index;
2265 
2266 private:
2267  template <typename Accessor>
2268  friend class TriaRawIterator;
2269  template <typename Accessor>
2270  friend class TriaIterator;
2271  template <typename Accessor>
2272  friend class TriaActiveIterator;
2273 };
2274 
2275 
2276 
2293 template <int spacedim>
2294 class TriaAccessor<0, 1, spacedim>
2295 {
2296 public:
2302  static const unsigned int space_dimension = spacedim;
2303 
2309  static const unsigned int dimension = 1;
2310 
2316  static const unsigned int structure_dimension = 0;
2317 
2321  using AccessorData = void;
2322 
2328  {
2340  right_vertex
2341  };
2342 
2355  const VertexKind vertex_kind,
2356  const unsigned int vertex_index);
2357 
2363  TriaAccessor(const Triangulation<1, spacedim> *tria = nullptr,
2364  const int = 0,
2365  const int = 0,
2366  const AccessorData * = nullptr);
2367 
2371  template <int structdim2, int dim2, int spacedim2>
2373 
2377  template <int structdim2, int dim2, int spacedim2>
2379 
2384  void
2385  copy_from(const TriaAccessor &);
2386 
2393  state();
2394 
2399  static int
2400  level();
2401 
2406  int
2407  index() const;
2408 
2414  get_triangulation() const;
2415 
2426  void
2427  operator++() const;
2428 
2433  void
2434  operator--() const;
2438  bool
2439  operator==(const TriaAccessor &) const;
2440 
2444  bool
2445  operator!=(const TriaAccessor &) const;
2446 
2455  bool
2456  operator<(const TriaAccessor &other) const;
2457 
2484  unsigned int
2485  vertex_index(const unsigned int i = 0) const;
2486 
2492  Point<spacedim> &
2493  vertex(const unsigned int i = 0) const;
2494 
2500  center() const;
2501 
2506  typename ::internal::TriangulationImplementation::
2507  Iterators<1, spacedim>::line_iterator static line(const unsigned int);
2508 
2515  static unsigned int
2516  line_index(const unsigned int i);
2517 
2521  static typename ::internal::TriangulationImplementation::
2522  Iterators<1, spacedim>::quad_iterator
2523  quad(const unsigned int i);
2524 
2531  static unsigned int
2532  quad_index(const unsigned int i);
2533 
2543  bool
2544  at_boundary() const;
2545 
2561  boundary_id() const;
2562 
2566  const Manifold<1, spacedim> &
2567  get_manifold() const;
2568 
2576  manifold_id() const;
2577 
2578 
2589  static unsigned char
2590  combined_face_orientation(const unsigned int face);
2591 
2595  static bool
2596  face_orientation(const unsigned int face);
2597 
2601  static bool
2602  face_flip(const unsigned int face);
2603 
2607  static bool
2608  face_rotation(const unsigned int face);
2609 
2613  static bool
2614  line_orientation(const unsigned int line);
2615 
2630  static bool
2631  has_children();
2632 
2637  static unsigned int
2638  n_children();
2639 
2644  static unsigned int
2645  n_active_descendants();
2646 
2651  static unsigned int
2652  number_of_children();
2653 
2654 
2658  static unsigned int
2659  max_refinement_depth();
2660 
2664  static unsigned int
2665  child_iterator_to_index(const TriaIterator<TriaAccessor<0, 1, spacedim>> &);
2666 
2671  child(const unsigned int);
2672 
2677  isotropic_child(const unsigned int);
2678 
2682  static RefinementCase<0>
2683  refinement_case();
2684 
2688  static int
2689  child_index(const unsigned int i);
2690 
2694  static int
2695  isotropic_child_index(const unsigned int i);
2728  void
2729  set_boundary_id(const types::boundary_id) const;
2730 
2737  void
2738  set_manifold_id(const types::manifold_id);
2739 
2751  void
2752  set_all_boundary_ids(const types::boundary_id) const;
2753 
2765  void
2766  set_all_manifold_ids(const types::manifold_id);
2774  bool
2775  used() const;
2776 
2781  reference_cell() const;
2782 
2786  unsigned int
2787  n_vertices() const;
2788 
2792  unsigned int
2793  n_lines() const;
2794 
2800  vertex_indices() const;
2801 
2807  line_indices() const;
2808 
2809 protected:
2814 
2820 
2824  unsigned int global_vertex_index;
2825 };
2826 
2827 
2828 
2844 template <int dim, int spacedim = dim>
2845 class CellAccessor : public TriaAccessor<dim, dim, spacedim>
2846 {
2847 public:
2852 
2857 
2868  CellAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
2869  const int level = -1,
2870  const int index = -1,
2871  const AccessorData * local_data = nullptr);
2872 
2876  CellAccessor(const TriaAccessor<dim, dim, spacedim> &cell_accessor);
2877 
2890  template <int structdim2, int dim2, int spacedim2>
2892 
2897  template <int structdim2, int dim2, int spacedim2>
2899 
2903  CellAccessor(const CellAccessor<dim, spacedim> &) = default;
2904 
2908  // NOLINTNEXTLINE OSX does not compile with noexcept
2910 
2914  ~CellAccessor() = default;
2915 
2926  operator=(const CellAccessor<dim, spacedim> &) = delete;
2927 
2931  // NOLINTNEXTLINE OSX does not compile with noexcept
2933  operator=(CellAccessor<dim, spacedim> &&) = default; // NOLINT
2934 
2951  child(const unsigned int i) const;
2952 
2956  boost::container::small_vector<TriaIterator<CellAccessor<dim, spacedim>>,
2958  child_iterators() const;
2959 
2963  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>
2964  face(const unsigned int i) const;
2965 
2970  unsigned int
2971  face_iterator_to_index(
2973 
2977  boost::container::small_vector<
2978  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>,
2980  face_iterators() const;
2981 
2991  unsigned int
2992  face_index(const unsigned int i) const;
2993 
3042  neighbor_child_on_subface(const unsigned int face_no,
3043  const unsigned int subface_no) const;
3044 
3070  neighbor(const unsigned int face_no) const;
3071 
3079  int
3080  neighbor_index(const unsigned int face_no) const;
3081 
3089  int
3090  neighbor_level(const unsigned int face_no) const;
3091 
3103  unsigned int
3104  neighbor_of_neighbor(const unsigned int face_no) const;
3105 
3116  bool
3117  neighbor_is_coarser(const unsigned int face_no) const;
3118 
3133  std::pair<unsigned int, unsigned int>
3134  neighbor_of_coarser_neighbor(const unsigned int neighbor) const;
3135 
3142  unsigned int
3143  neighbor_face_no(const unsigned int neighbor) const;
3144 
3148  static bool
3149  is_level_cell();
3150 
3164  bool
3165  has_periodic_neighbor(const unsigned int i) const;
3166 
3184  periodic_neighbor(const unsigned int i) const;
3185 
3194  neighbor_or_periodic_neighbor(const unsigned int i) const;
3195 
3211  periodic_neighbor_child_on_subface(const unsigned int face_no,
3212  const unsigned int subface_no) const;
3213 
3224  std::pair<unsigned int, unsigned int>
3225  periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const;
3226 
3232  int
3233  periodic_neighbor_index(const unsigned int i) const;
3234 
3240  int
3241  periodic_neighbor_level(const unsigned int i) const;
3242 
3257  unsigned int
3258  periodic_neighbor_of_periodic_neighbor(const unsigned int i) const;
3259 
3265  unsigned int
3266  periodic_neighbor_face_no(const unsigned int i) const;
3267 
3274  bool
3275  periodic_neighbor_is_coarser(const unsigned int i) const;
3276 
3293  bool
3294  at_boundary(const unsigned int i) const;
3295 
3304  bool
3305  at_boundary() const;
3306 
3314  bool
3315  has_boundary_lines() const;
3342  refine_flag_set() const;
3343 
3361  void
3362  set_refine_flag(const RefinementCase<dim> ref_case =
3364 
3368  void
3369  clear_refine_flag() const;
3370 
3378  bool
3379  flag_for_face_refinement(
3380  const unsigned int face_no,
3381  const RefinementCase<dim - 1> &face_refinement_case =
3383 
3389  bool
3390  flag_for_line_refinement(const unsigned int line_no) const;
3391 
3401  subface_case(const unsigned int face_no) const;
3402 
3406  bool
3407  coarsen_flag_set() const;
3408 
3413  void
3414  set_coarsen_flag() const;
3415 
3419  void
3420  clear_coarsen_flag() const;
3444  material_id() const;
3445 
3457  void
3458  set_material_id(const types::material_id new_material_id) const;
3459 
3468  void
3469  recursively_set_material_id(const types::material_id new_material_id) const;
3496  subdomain_id() const;
3497 
3513  void
3514  set_subdomain_id(const types::subdomain_id new_subdomain_id) const;
3515 
3525  level_subdomain_id() const;
3526 
3531  void
3532  set_level_subdomain_id(
3533  const types::subdomain_id new_level_subdomain_id) const;
3534 
3535 
3551  void
3552  recursively_set_subdomain_id(
3553  const types::subdomain_id new_subdomain_id) const;
3576  global_active_cell_index() const;
3577 
3584  global_level_cell_index() const;
3585 
3599  bool
3600  direction_flag() const;
3601 
3627  unsigned int
3628  active_cell_index() const;
3629 
3637  int
3638  parent_index() const;
3639 
3646  parent() const;
3647 
3667  bool
3668  is_active() const;
3669 
3689  bool
3690  is_locally_owned() const;
3691 
3696  bool
3697  is_locally_owned_on_level() const;
3698 
3732  bool
3733  is_ghost() const;
3734 
3740  bool
3741  is_ghost_on_level() const;
3742 
3769  bool
3770  is_artificial() const;
3771 
3778  bool
3779  is_artificial_on_level() const;
3780 
3794  bool
3795  point_inside(const Point<spacedim> &p) const;
3796 
3805  void
3806  set_neighbor(const unsigned int i,
3807  const TriaIterator<CellAccessor<dim, spacedim>> &pointer) const;
3808 
3822  CellId
3823  id() const;
3824 
3826 
3830  double
3831  diameter(const Mapping<dim, spacedim> &mapping) const;
3832 
3841  DeclException0(ExcRefineCellNotActive);
3845  DeclException0(ExcCellFlaggedForRefinement);
3849  DeclException0(ExcCellFlaggedForCoarsening);
3850 
3851 protected:
3867  unsigned int
3868  neighbor_of_neighbor_internal(const unsigned int neighbor) const;
3869 
3875  template <int dim_, int spacedim_>
3876  bool
3877  point_inside_codim(const Point<spacedim_> &p) const;
3878 
3879 
3880 
3881 private:
3886  void
3887  set_active_cell_index(const unsigned int active_cell_index) const;
3888 
3892  void
3893  set_global_active_cell_index(const types::global_cell_index index) const;
3894 
3898  void
3899  set_global_level_cell_index(const types::global_cell_index index) const;
3900 
3904  void
3905  set_parent(const unsigned int parent_index);
3906 
3913  void
3914  set_direction_flag(const bool new_direction_flag) const;
3915 
3916  template <int, int>
3917  friend class Triangulation;
3918 
3919  template <int, int>
3921 
3922  friend struct ::internal::TriangulationImplementation::Implementation;
3923  friend struct ::internal::TriangulationImplementation::
3924  ImplementationMixedMesh;
3925 };
3926 
3927 
3928 
3929 /* ----- declaration of explicit specializations and general templates ----- */
3930 
3931 
3932 template <int structdim, int dim, int spacedim>
3933 template <typename OtherAccessor>
3935  const OtherAccessor &)
3936 {
3937  Assert(false,
3938  ExcMessage("You are attempting an illegal conversion between "
3939  "iterator/accessor types. The constructor you call "
3940  "only exists to make certain template constructs "
3941  "easier to write as dimension independent code but "
3942  "the conversion is not valid in the current context."));
3943 }
3944 
3945 
3946 
3947 template <int structdim, int dim, int spacedim>
3948 template <int structdim2, int dim2, int spacedim2>
3951 {
3952  Assert(false,
3953  ExcMessage("You are attempting an illegal conversion between "
3954  "iterator/accessor types. The constructor you call "
3955  "only exists to make certain template constructs "
3956  "easier to write as dimension independent code but "
3957  "the conversion is not valid in the current context."));
3958 }
3959 
3960 
3961 
3962 template <int dim, int spacedim>
3963 template <int structdim2, int dim2, int spacedim2>
3966 {
3967  Assert(false,
3968  ExcMessage("You are attempting an illegal conversion between "
3969  "iterator/accessor types. The constructor you call "
3970  "only exists to make certain template constructs "
3971  "easier to write as dimension independent code but "
3972  "the conversion is not valid in the current context."));
3973 }
3974 
3975 
3976 
3977 template <int structdim, int dim, int spacedim>
3978 template <int structdim2, int dim2, int spacedim2>
3981 {
3982  Assert(false,
3983  ExcMessage("You are attempting an illegal conversion between "
3984  "iterator/accessor types. The constructor you call "
3985  "only exists to make certain template constructs "
3986  "easier to write as dimension independent code but "
3987  "the conversion is not valid in the current context."));
3988 }
3989 
3990 
3991 
3992 template <int dim, int spacedim>
3993 template <int structdim2, int dim2, int spacedim2>
3996 {
3997  Assert(false,
3998  ExcMessage("You are attempting an illegal conversion between "
3999  "iterator/accessor types. The constructor you call "
4000  "only exists to make certain template constructs "
4001  "easier to write as dimension independent code but "
4002  "the conversion is not valid in the current context."));
4003 }
4004 
4005 
4006 #ifndef DOXYGEN
4007 
4008 template <>
4009 bool
4011 template <>
4012 bool
4014 template <>
4015 bool
4017 template <>
4018 bool
4020 template <>
4021 bool
4023 template <>
4024 bool
4026 // -------------------------------------------------------------------
4027 
4028 template <>
4029 void
4031 
4032 #endif // DOXYGEN
4033 
4035 
4036 // include more templates in debug and optimized mode
4037 #include "tria_accessor.templates.h"
4038 
4039 #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:82
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< dim, spacedim > & tria
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:414
static ::ExceptionBase & ExcMessage(std::string arg1)
unsigned int subdomain_id
Definition: types.h:43
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:509
const Triangulation< dim, spacedim > * tria
static ::ExceptionBase & ExcSetOnlyEvenChildren(int arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1461
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:487
unsigned int vertex_indices[2]
Definition: grid_tools.cc:1256
#define DeclException0(Exception0)
Definition: exceptions.h:464
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
unsigned int level
Definition: grid_out.cc:4589
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:185
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()
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
bool point_inside(const Point< spacedim > &p) const
#define DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
Definition: config.h:452
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
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:160
unsigned int boundary_id
Definition: types.h:129
typename ::internal::TriaAccessorImplementation::PresentLevelType< structdim, dim >::type present_level
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcCantCompareIterators()