Reference documentation for deal.II version Git 9e557027ad 2021-09-25 18:07:42 +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;
913  bool
914  face_orientation(const unsigned int face) const;
915 
925  bool
926  face_flip(const unsigned int face) const;
927 
937  bool
938  face_rotation(const unsigned int face) const;
939 
952  bool
953  line_orientation(const unsigned int line) const;
968  bool
969  has_children() const;
970 
975  unsigned int
976  n_children() const;
977 
982  unsigned int
983  number_of_children() const;
984 
998  unsigned int
999  n_active_descendants() const;
1000 
1014  unsigned int
1015  max_refinement_depth() const;
1016 
1021  child(const unsigned int i) const;
1022 
1027  unsigned int
1028  child_iterator_to_index(
1030 
1040  isotropic_child(const unsigned int i) const;
1041 
1046  refinement_case() const;
1047 
1053  int
1054  child_index(const unsigned int i) const;
1055 
1061  int
1062  isotropic_child_index(const unsigned int i) const;
1085  boundary_id() const;
1086 
1116  void
1117  set_boundary_id(const types::boundary_id) const;
1118 
1149  void
1150  set_all_boundary_ids(const types::boundary_id) const;
1151 
1159  bool
1160  at_boundary() const;
1161 
1171  const Manifold<dim, spacedim> &
1172  get_manifold() const;
1173 
1195  manifold_id() const;
1196 
1214  void
1215  set_manifold_id(const types::manifold_id) const;
1216 
1230  void
1231  set_all_manifold_ids(const types::manifold_id) const;
1232 
1249  bool
1250  user_flag_set() const;
1251 
1257  void
1258  set_user_flag() const;
1259 
1265  void
1266  clear_user_flag() const;
1267 
1273  void
1274  recursively_set_user_flag() const;
1275 
1281  void
1282  recursively_clear_user_flag() const;
1283 
1289  void
1290  clear_user_data() const;
1291 
1303  void
1304  set_user_pointer(void *p) const;
1305 
1311  void
1312  clear_user_pointer() const;
1313 
1329  void *
1330  user_pointer() const;
1331 
1353  void
1354  recursively_set_user_pointer(void *p) const;
1355 
1362  void
1363  recursively_clear_user_pointer() const;
1364 
1374  void
1375  set_user_index(const unsigned int p) const;
1376 
1382  void
1383  clear_user_index() const;
1384 
1396  unsigned int
1397  user_index() const;
1398 
1416  void
1417  recursively_set_user_index(const unsigned int p) const;
1418 
1427  void
1428  recursively_clear_user_index() const;
1463  double
1464  diameter() const;
1465 
1492  std::pair<Point<spacedim>, double>
1493  enclosing_ball() const;
1494 
1504  bounding_box() const;
1505 
1515  double
1516  extent_in_direction(const unsigned int axis) const;
1517 
1521  double
1522  minimum_vertex_distance() const;
1523 
1538  intermediate_point(const Point<structdim> &coordinates) const;
1539 
1563  real_to_unit_cell_affine_approximation(const Point<spacedim> &point) const;
1564 
1600  center(const bool respect_manifold = false,
1601  const bool interpolate_from_surrounding = false) const;
1602 
1621  barycenter() const;
1622 
1648  double
1649  measure() const;
1650 
1665  bool
1666  is_translation_of(
1668 
1673  reference_cell() const;
1674 
1678  unsigned int
1679  n_vertices() const;
1680 
1684  unsigned int
1685  n_lines() const;
1686 
1692  unsigned int
1693  n_faces() const;
1694 
1700  vertex_indices() const;
1701 
1707  line_indices() const;
1708 
1716  face_indices() const;
1717 
1722 private:
1727  void
1728  set_boundary_id_internal(const types::boundary_id id) const;
1729 
1737  void
1738  set_bounding_object_indices(
1739  const std::initializer_list<int> &new_indices) const;
1740 
1744  void
1745  set_bounding_object_indices(
1746  const std::initializer_list<unsigned int> &new_indices) const;
1747 
1755  void
1756  set_line_orientation(const unsigned int line, const bool orientation) const;
1757 
1768  void
1769  set_face_orientation(const unsigned int face, const bool orientation) const;
1770 
1777  void
1778  set_face_flip(const unsigned int face, const bool flip) const;
1779 
1786  void
1787  set_face_rotation(const unsigned int face, const bool rotation) const;
1788 
1792  void
1793  set_used_flag() const;
1794 
1798  void
1799  clear_used_flag() const;
1800 
1809  void
1810  set_refinement_case(const RefinementCase<structdim> &ref_case) const;
1811 
1819  void
1820  clear_refinement_case() const;
1821 
1828  void
1829  set_children(const unsigned int i, const int index) const;
1830 
1835  void
1836  clear_children() const;
1837 
1838 private:
1839  template <int, int>
1840  friend class Triangulation;
1841 
1842  friend struct ::internal::TriangulationImplementation::Implementation;
1843  friend struct ::internal::TriangulationImplementation::
1844  ImplementationMixedMesh;
1845  friend struct ::internal::TriaAccessorImplementation::Implementation;
1846 };
1847 
1848 
1849 
1868 template <int dim, int spacedim>
1869 class TriaAccessor<0, dim, spacedim>
1870 {
1871 public:
1877  static const unsigned int space_dimension = spacedim;
1878 
1884  static const unsigned int dimension = dim;
1885 
1891  static const unsigned int structure_dimension = 0;
1892 
1896  using AccessorData = void;
1897 
1903  const unsigned int vertex_index);
1904 
1911  const int level = 0,
1912  const int index = 0,
1913  const AccessorData * = nullptr);
1914 
1918  template <int structdim2, int dim2, int spacedim2>
1920 
1924  template <int structdim2, int dim2, int spacedim2>
1926 
1931  state() const;
1932 
1937  static int
1938  level();
1939 
1944  int
1945  index() const;
1946 
1952  get_triangulation() const;
1953 
1963  void
1964  operator++();
1965 
1969  void
1970  operator--();
1974  bool
1975  operator==(const TriaAccessor &) const;
1976 
1980  bool
1981  operator!=(const TriaAccessor &) const;
1982 
2010  unsigned int
2011  vertex_index(const unsigned int i = 0) const;
2012 
2018  Point<spacedim> &
2019  vertex(const unsigned int i = 0) const;
2020 
2025  typename ::internal::TriangulationImplementation::
2026  Iterators<dim, spacedim>::line_iterator static line(const unsigned int);
2027 
2031  static unsigned int
2032  line_index(const unsigned int i);
2033 
2037  static typename ::internal::TriangulationImplementation::
2038  Iterators<dim, spacedim>::quad_iterator
2039  quad(const unsigned int i);
2040 
2044  static unsigned int
2045  quad_index(const unsigned int i);
2046 
2062  double
2063  diameter() const;
2064 
2072  double
2073  extent_in_direction(const unsigned int axis) const;
2074 
2083  center(const bool respect_manifold = false,
2084  const bool interpolate_from_surrounding = false) const;
2085 
2094  double
2095  measure() const;
2110  static bool
2111  face_orientation(const unsigned int face);
2112 
2116  static bool
2117  face_flip(const unsigned int face);
2118 
2122  static bool
2123  face_rotation(const unsigned int face);
2124 
2128  static bool
2129  line_orientation(const unsigned int line);
2130 
2145  static bool
2146  has_children();
2147 
2152  static unsigned int
2153  n_children();
2154 
2159  static unsigned int
2160  n_active_descendants();
2161 
2166  static unsigned int
2167  number_of_children();
2168 
2169 
2173  static unsigned int
2174  max_refinement_depth();
2175 
2179  static unsigned int
2180  child_iterator_to_index(const TriaIterator<TriaAccessor<0, dim, spacedim>> &);
2181 
2186  child(const unsigned int);
2187 
2192  isotropic_child(const unsigned int);
2193 
2197  static RefinementCase<0>
2198  refinement_case();
2199 
2203  static int
2204  child_index(const unsigned int i);
2205 
2209  static int
2210  isotropic_child_index(const unsigned int i);
2218  bool
2219  used() const;
2220 
2221 protected:
2229  void
2230  copy_from(const TriaAccessor &);
2231 
2240  bool
2241  operator<(const TriaAccessor &other) const;
2242 
2247 
2251  unsigned int global_vertex_index;
2252 
2253 private:
2254  template <typename Accessor>
2255  friend class TriaRawIterator;
2256  template <typename Accessor>
2257  friend class TriaIterator;
2258  template <typename Accessor>
2259  friend class TriaActiveIterator;
2260 };
2261 
2262 
2263 
2280 template <int spacedim>
2281 class TriaAccessor<0, 1, spacedim>
2282 {
2283 public:
2289  static const unsigned int space_dimension = spacedim;
2290 
2296  static const unsigned int dimension = 1;
2297 
2303  static const unsigned int structure_dimension = 0;
2304 
2308  using AccessorData = void;
2309 
2315  {
2327  right_vertex
2328  };
2329 
2342  const VertexKind vertex_kind,
2343  const unsigned int vertex_index);
2344 
2350  TriaAccessor(const Triangulation<1, spacedim> *tria = nullptr,
2351  const int = 0,
2352  const int = 0,
2353  const AccessorData * = nullptr);
2354 
2358  template <int structdim2, int dim2, int spacedim2>
2360 
2364  template <int structdim2, int dim2, int spacedim2>
2366 
2371  void
2372  copy_from(const TriaAccessor &);
2373 
2380  state();
2381 
2386  static int
2387  level();
2388 
2393  int
2394  index() const;
2395 
2401  get_triangulation() const;
2402 
2413  void
2414  operator++() const;
2415 
2420  void
2421  operator--() const;
2425  bool
2426  operator==(const TriaAccessor &) const;
2427 
2431  bool
2432  operator!=(const TriaAccessor &) const;
2433 
2442  bool
2443  operator<(const TriaAccessor &other) const;
2444 
2471  unsigned int
2472  vertex_index(const unsigned int i = 0) const;
2473 
2479  Point<spacedim> &
2480  vertex(const unsigned int i = 0) const;
2481 
2487  center() const;
2488 
2493  typename ::internal::TriangulationImplementation::
2494  Iterators<1, spacedim>::line_iterator static line(const unsigned int);
2495 
2502  static unsigned int
2503  line_index(const unsigned int i);
2504 
2508  static typename ::internal::TriangulationImplementation::
2509  Iterators<1, spacedim>::quad_iterator
2510  quad(const unsigned int i);
2511 
2518  static unsigned int
2519  quad_index(const unsigned int i);
2520 
2530  bool
2531  at_boundary() const;
2532 
2548  boundary_id() const;
2549 
2553  const Manifold<1, spacedim> &
2554  get_manifold() const;
2555 
2563  manifold_id() const;
2564 
2565 
2576  static bool
2577  face_orientation(const unsigned int face);
2578 
2582  static bool
2583  face_flip(const unsigned int face);
2584 
2588  static bool
2589  face_rotation(const unsigned int face);
2590 
2594  static bool
2595  line_orientation(const unsigned int line);
2596 
2611  static bool
2612  has_children();
2613 
2618  static unsigned int
2619  n_children();
2620 
2625  static unsigned int
2626  n_active_descendants();
2627 
2632  static unsigned int
2633  number_of_children();
2634 
2635 
2639  static unsigned int
2640  max_refinement_depth();
2641 
2645  static unsigned int
2646  child_iterator_to_index(const TriaIterator<TriaAccessor<0, 1, spacedim>> &);
2647 
2652  child(const unsigned int);
2653 
2658  isotropic_child(const unsigned int);
2659 
2663  static RefinementCase<0>
2664  refinement_case();
2665 
2669  static int
2670  child_index(const unsigned int i);
2671 
2675  static int
2676  isotropic_child_index(const unsigned int i);
2709  void
2710  set_boundary_id(const types::boundary_id) const;
2711 
2718  void
2719  set_manifold_id(const types::manifold_id);
2720 
2732  void
2733  set_all_boundary_ids(const types::boundary_id) const;
2734 
2746  void
2747  set_all_manifold_ids(const types::manifold_id);
2755  bool
2756  used() const;
2757 
2762  reference_cell() const;
2763 
2767  unsigned int
2768  n_vertices() const;
2769 
2773  unsigned int
2774  n_lines() const;
2775 
2781  vertex_indices() const;
2782 
2788  line_indices() const;
2789 
2790 protected:
2795 
2801 
2805  unsigned int global_vertex_index;
2806 };
2807 
2808 
2809 
2825 template <int dim, int spacedim = dim>
2826 class CellAccessor : public TriaAccessor<dim, dim, spacedim>
2827 {
2828 public:
2833 
2838 
2849  CellAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
2850  const int level = -1,
2851  const int index = -1,
2852  const AccessorData * local_data = nullptr);
2853 
2857  CellAccessor(const TriaAccessor<dim, dim, spacedim> &cell_accessor);
2858 
2871  template <int structdim2, int dim2, int spacedim2>
2873 
2878  template <int structdim2, int dim2, int spacedim2>
2880 
2884  CellAccessor(const CellAccessor<dim, spacedim> &) = default;
2885 
2889  // NOLINTNEXTLINE OSX does not compile with noexcept
2891 
2895  ~CellAccessor() = default;
2896 
2907  operator=(const CellAccessor<dim, spacedim> &) = delete;
2908 
2912  // NOLINTNEXTLINE OSX does not compile with noexcept
2914  operator=(CellAccessor<dim, spacedim> &&) = default; // NOLINT
2915 
2932  child(const unsigned int i) const;
2933 
2937  boost::container::small_vector<TriaIterator<CellAccessor<dim, spacedim>>,
2939  child_iterators() const;
2940 
2944  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>
2945  face(const unsigned int i) const;
2946 
2951  unsigned int
2952  face_iterator_to_index(
2954 
2958  boost::container::small_vector<
2959  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>,
2961  face_iterators() const;
2962 
2972  unsigned int
2973  face_index(const unsigned int i) const;
2974 
3023  neighbor_child_on_subface(const unsigned int face_no,
3024  const unsigned int subface_no) const;
3025 
3051  neighbor(const unsigned int face_no) const;
3052 
3060  int
3061  neighbor_index(const unsigned int face_no) const;
3062 
3070  int
3071  neighbor_level(const unsigned int face_no) const;
3072 
3084  unsigned int
3085  neighbor_of_neighbor(const unsigned int face_no) const;
3086 
3097  bool
3098  neighbor_is_coarser(const unsigned int face_no) const;
3099 
3114  std::pair<unsigned int, unsigned int>
3115  neighbor_of_coarser_neighbor(const unsigned int neighbor) const;
3116 
3123  unsigned int
3124  neighbor_face_no(const unsigned int neighbor) const;
3125 
3129  static bool
3130  is_level_cell();
3131 
3145  bool
3146  has_periodic_neighbor(const unsigned int i) const;
3147 
3165  periodic_neighbor(const unsigned int i) const;
3166 
3175  neighbor_or_periodic_neighbor(const unsigned int i) const;
3176 
3192  periodic_neighbor_child_on_subface(const unsigned int face_no,
3193  const unsigned int subface_no) const;
3194 
3205  std::pair<unsigned int, unsigned int>
3206  periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const;
3207 
3213  int
3214  periodic_neighbor_index(const unsigned int i) const;
3215 
3221  int
3222  periodic_neighbor_level(const unsigned int i) const;
3223 
3238  unsigned int
3239  periodic_neighbor_of_periodic_neighbor(const unsigned int i) const;
3240 
3246  unsigned int
3247  periodic_neighbor_face_no(const unsigned int i) const;
3248 
3255  bool
3256  periodic_neighbor_is_coarser(const unsigned int i) const;
3257 
3274  bool
3275  at_boundary(const unsigned int i) const;
3276 
3285  bool
3286  at_boundary() const;
3287 
3295  bool
3296  has_boundary_lines() const;
3323  refine_flag_set() const;
3324 
3342  void
3343  set_refine_flag(const RefinementCase<dim> ref_case =
3345 
3349  void
3350  clear_refine_flag() const;
3351 
3359  bool
3360  flag_for_face_refinement(
3361  const unsigned int face_no,
3362  const RefinementCase<dim - 1> &face_refinement_case =
3364 
3370  bool
3371  flag_for_line_refinement(const unsigned int line_no) const;
3372 
3382  subface_case(const unsigned int face_no) const;
3383 
3387  bool
3388  coarsen_flag_set() const;
3389 
3394  void
3395  set_coarsen_flag() const;
3396 
3400  void
3401  clear_coarsen_flag() const;
3425  material_id() const;
3426 
3438  void
3439  set_material_id(const types::material_id new_material_id) const;
3440 
3449  void
3450  recursively_set_material_id(const types::material_id new_material_id) const;
3477  subdomain_id() const;
3478 
3494  void
3495  set_subdomain_id(const types::subdomain_id new_subdomain_id) const;
3496 
3506  level_subdomain_id() const;
3507 
3512  void
3513  set_level_subdomain_id(
3514  const types::subdomain_id new_level_subdomain_id) const;
3515 
3516 
3532  void
3533  recursively_set_subdomain_id(
3534  const types::subdomain_id new_subdomain_id) const;
3557  global_active_cell_index() const;
3558 
3565  global_level_cell_index() const;
3566 
3580  bool
3581  direction_flag() const;
3582 
3608  unsigned int
3609  active_cell_index() const;
3610 
3618  int
3619  parent_index() const;
3620 
3627  parent() const;
3628 
3648  bool
3649  is_active() const;
3650 
3670  bool
3671  is_locally_owned() const;
3672 
3677  bool
3678  is_locally_owned_on_level() const;
3679 
3713  bool
3714  is_ghost() const;
3715 
3721  bool
3722  is_ghost_on_level() const;
3723 
3750  bool
3751  is_artificial() const;
3752 
3759  bool
3760  is_artificial_on_level() const;
3761 
3775  bool
3776  point_inside(const Point<spacedim> &p) const;
3777 
3786  void
3787  set_neighbor(const unsigned int i,
3788  const TriaIterator<CellAccessor<dim, spacedim>> &pointer) const;
3789 
3803  CellId
3804  id() const;
3805 
3807 
3811  double
3812  diameter(const Mapping<dim, spacedim> &mapping) const;
3813 
3822  DeclException0(ExcRefineCellNotActive);
3826  DeclException0(ExcCellFlaggedForRefinement);
3830  DeclException0(ExcCellFlaggedForCoarsening);
3831 
3832 protected:
3848  unsigned int
3849  neighbor_of_neighbor_internal(const unsigned int neighbor) const;
3850 
3856  template <int dim_, int spacedim_>
3857  bool
3858  point_inside_codim(const Point<spacedim_> &p) const;
3859 
3860 
3861 
3862 private:
3867  void
3868  set_active_cell_index(const unsigned int active_cell_index) const;
3869 
3873  void
3874  set_global_active_cell_index(const types::global_cell_index index) const;
3875 
3879  void
3880  set_global_level_cell_index(const types::global_cell_index index) const;
3881 
3885  void
3886  set_parent(const unsigned int parent_index);
3887 
3894  void
3895  set_direction_flag(const bool new_direction_flag) const;
3896 
3897  template <int, int>
3898  friend class Triangulation;
3899 
3900  template <int, int>
3902 
3903  friend struct ::internal::TriangulationImplementation::Implementation;
3904  friend struct ::internal::TriangulationImplementation::
3905  ImplementationMixedMesh;
3906 };
3907 
3908 
3909 
3910 /* ----- declaration of explicit specializations and general templates ----- */
3911 
3912 
3913 template <int structdim, int dim, int spacedim>
3914 template <typename OtherAccessor>
3916  const OtherAccessor &)
3917 {
3918  Assert(false,
3919  ExcMessage("You are attempting an illegal conversion between "
3920  "iterator/accessor types. The constructor you call "
3921  "only exists to make certain template constructs "
3922  "easier to write as dimension independent code but "
3923  "the conversion is not valid in the current context."));
3924 }
3925 
3926 
3927 
3928 template <int structdim, int dim, int spacedim>
3929 template <int structdim2, int dim2, int spacedim2>
3932 {
3933  Assert(false,
3934  ExcMessage("You are attempting an illegal conversion between "
3935  "iterator/accessor types. The constructor you call "
3936  "only exists to make certain template constructs "
3937  "easier to write as dimension independent code but "
3938  "the conversion is not valid in the current context."));
3939 }
3940 
3941 
3942 
3943 template <int dim, int spacedim>
3944 template <int structdim2, int dim2, int spacedim2>
3947 {
3948  Assert(false,
3949  ExcMessage("You are attempting an illegal conversion between "
3950  "iterator/accessor types. The constructor you call "
3951  "only exists to make certain template constructs "
3952  "easier to write as dimension independent code but "
3953  "the conversion is not valid in the current context."));
3954 }
3955 
3956 
3957 
3958 template <int structdim, int dim, int spacedim>
3959 template <int structdim2, int dim2, int spacedim2>
3962 {
3963  Assert(false,
3964  ExcMessage("You are attempting an illegal conversion between "
3965  "iterator/accessor types. The constructor you call "
3966  "only exists to make certain template constructs "
3967  "easier to write as dimension independent code but "
3968  "the conversion is not valid in the current context."));
3969 }
3970 
3971 
3972 
3973 template <int dim, int spacedim>
3974 template <int structdim2, int dim2, int spacedim2>
3977 {
3978  Assert(false,
3979  ExcMessage("You are attempting an illegal conversion between "
3980  "iterator/accessor types. The constructor you call "
3981  "only exists to make certain template constructs "
3982  "easier to write as dimension independent code but "
3983  "the conversion is not valid in the current context."));
3984 }
3985 
3986 
3987 #ifndef DOXYGEN
3988 
3989 template <>
3990 bool
3992 template <>
3993 bool
3995 template <>
3996 bool
3998 template <>
3999 bool
4001 template <>
4002 bool
4004 template <>
4005 bool
4007 // -------------------------------------------------------------------
4008 
4009 template <>
4010 void
4012 
4013 #endif // DOXYGEN
4014 
4016 
4017 // include more templates in debug and optimized mode
4018 #include "tria_accessor.templates.h"
4019 
4020 #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< 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:1254
#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()