Reference documentation for deal.II version Git a5ed68a04a 2019-09-22 06:50:58 -0600
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
tria_accessor.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2019 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 
22 #include <deal.II/base/bounding_box.h>
23 #include <deal.II/base/exceptions.h>
24 #include <deal.II/base/geometry_info.h>
25 #include <deal.II/base/point.h>
26 
27 #include <deal.II/grid/cell_id.h>
28 #include <deal.II/grid/tria_iterator_base.h>
29 #include <deal.II/grid/tria_iterator_selector.h>
30 
31 #include <utility>
32 
33 
34 DEAL_II_NAMESPACE_OPEN
35 
36 // Forward declarations
37 #ifndef DOXYGEN
38 template <int dim, int spacedim>
39 class Triangulation;
40 template <typename Accessor>
41 class TriaRawIterator;
42 template <typename Accessor>
43 class TriaIterator;
44 template <typename Accessor>
45 class TriaActiveIterator;
46 
47 template <int dim, int spacedim>
48 class Manifold;
49 #endif
50 
51 namespace internal
52 {
53  namespace TriangulationImplementation
54  {
55  template <int dim>
56  class TriaObject;
57  template <typename G>
58  class TriaObjects;
59  struct Implementation;
60  } // namespace TriangulationImplementation
61 
62  namespace TriaAccessorImplementation
63  {
64  struct Implementation;
65 
71  template <int structdim, int dim>
73  {
74  struct type
75  {
79  type() = default;
80 
84  type(const int level)
85  {
86  Assert(level == 0, ExcInternalError());
87  (void)level; // removes -Wunused-parameter warning in optimized mode
88  }
89 
93  operator int() const
94  {
95  return 0;
96  }
97 
98  void
99  operator++() const
100  {
101  Assert(false, ExcInternalError());
102  }
103 
104  void
105  operator--() const
106  {
107  Assert(false, ExcInternalError());
108  }
109  };
110  };
111 
112 
118  template <int dim>
119  struct PresentLevelType<dim, dim>
120  {
121  using type = int;
122  };
123 
124  } // namespace TriaAccessorImplementation
125 } // namespace internal
126 template <int structdim, int dim, int spacedim>
128 template <int dim, int spacedim>
129 class TriaAccessor<0, dim, spacedim>;
130 template <int spacedim>
131 class TriaAccessor<0, 1, spacedim>;
132 
133 // note: the file tria_accessor.templates.h is included at the end of
134 // this file. this includes a lot of templates. originally, this was
135 // only done in debug mode, but led to cyclic reduction problems and
136 // so is now on by default.
137 
138 
143 {
148  "The operation you are attempting can only be performed for "
149  "(cell, face, or edge) iterators that point to valid "
150  "objects. These objects need not necessarily be active, "
151  "i.e., have no children, but they need to be part of a "
152  "triangulation. (The objects pointed to by an iterator "
153  "may -- after coarsening -- also be objects that used "
154  "to be part of a triangulation, but are now no longer "
155  "used. Their memory location may have been retained "
156  "for re-use upon the next mesh refinement, but is "
157  "currently unused.)");
168  "The operation you are attempting can only be performed for "
169  "(cell, face, or edge) iterators that point to 'active' "
170  "objects. 'Active' objects are those that do not have "
171  "children (in the case of cells), or that are part of "
172  "an active cell (in the case of faces or edges). However, "
173  "the object on which you are trying the current "
174  "operation is not 'active' in this sense.");
181  "The operation you are attempting can only be performed for "
182  "(cell, face, or edge) iterators that have children, "
183  "but the object on which you are trying the current "
184  "operation does not have any.");
192  "The operation you are attempting can only be performed for "
193  "(cell, face, or edge) iterators that have a parent object, "
194  "but the object on which you are trying the current "
195  "operation does not have one -- i.e., it is on the "
196  "coarsest level of the triangulation.");
201  int,
202  << "You can only set the child index if the cell does not "
203  << "currently have children registered; or you can clear it. "
204  << "The given index was " << arg1
205  << " (-1 means: clear children).");
209  template <typename AccessorType>
211  AccessorType,
212  << "You tried to dereference an iterator for which this "
213  << "is not possible. More information on this iterator: "
214  << "index=" << arg1.index() << ", state="
215  << (arg1.state() == IteratorState::valid ?
216  "valid" :
217  (arg1.state() == IteratorState::past_the_end ?
218  "past_the_end" :
219  "invalid")));
224  "Iterators can only be compared if they point to the same "
225  "triangulation, or if neither of them are associated "
226  "with a triangulation.");
227  // TODO: Write documentation!
232  // TODO: Write documentation!
252  // TODO: Write documentation!
258  int,
259  << "You can only set the child index of an even numbered child."
260  << "The number of the child given was " << arg1 << ".");
261 } // namespace TriaAccessorExceptions
262 
263 
290 template <int structdim, int dim, int spacedim = dim>
292 {
293 public:
299  static const unsigned int space_dimension = spacedim;
300 
306  static const unsigned int dimension = dim;
307 
313  static const unsigned int structure_dimension = structdim;
314 
324  void
325  operator=(const TriaAccessorBase *) = delete;
326 
327 protected:
333  using AccessorData = void;
334 
338  TriaAccessorBase(const Triangulation<dim, spacedim> *parent = nullptr,
339  const int level = -1,
340  const int index = -1,
341  const AccessorData * = nullptr);
342 
347 
355  void
356  copy_from(const TriaAccessorBase &);
357 
362  operator=(const TriaAccessorBase &);
363 
374  bool
375  operator<(const TriaAccessorBase &other) const;
376 
377 protected:
381  bool
382  operator==(const TriaAccessorBase &) const;
383 
387  bool
388  operator!=(const TriaAccessorBase &) const;
389 
403  void
404  operator++();
405 
413  void
414  operator--();
424  objects() const;
425 
426 public:
432  using LocalData = void *;
433 
457  int
458  level() const;
459 
486  int
487  index() const;
488 
494  state() const;
495 
501  get_triangulation() const;
502 
506 protected:
511  typename ::internal::TriaAccessorImplementation::
512  PresentLevelType<structdim, dim>::type present_level;
513 
519 
524 
525 private:
526  template <typename Accessor>
527  friend class TriaRawIterator;
528  template <typename Accessor>
529  friend class TriaIterator;
530  template <typename Accessor>
531  friend class TriaActiveIterator;
532 };
533 
534 
535 
557 template <int structdim, int dim, int spacedim = dim>
558 class InvalidAccessor : public TriaAccessorBase<structdim, dim, spacedim>
559 {
560 public:
564  using AccessorData =
566 
574  InvalidAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
575  const int level = -1,
576  const int index = -1,
577  const AccessorData * local_data = nullptr);
578 
587 
592  template <typename OtherAccessor>
593  InvalidAccessor(const OtherAccessor &);
594 
598  void
599  copy_from(const InvalidAccessor &);
600 
604  bool
605  operator==(const InvalidAccessor &) const;
606  bool
607  operator!=(const InvalidAccessor &) const;
608 
612  void
613  operator++() const;
614  void
615  operator--() const;
616 
621  bool
622  used() const;
623 
628  bool
629  has_children() const;
630 
635  manifold_id() const;
636 
640  unsigned int
641  user_index() const;
642 
646  void
647  set_user_index(const unsigned int p) const;
648 
652  void
653  set_manifold_id(const types::manifold_id) const;
654 
659  vertex(const unsigned int i) const;
660 
665  typename ::internal::TriangulationImplementation::
666  Iterators<dim, spacedim>::line_iterator
667  line(const unsigned int i) const;
668 
673  typename ::internal::TriangulationImplementation::
674  Iterators<dim, spacedim>::quad_iterator
675  quad(const unsigned int i) const;
676 };
677 
678 
679 
698 template <int structdim, int dim, int spacedim>
699 class TriaAccessor : public TriaAccessorBase<structdim, dim, spacedim>
700 {
701 public:
705  using AccessorData =
707 
711  TriaAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
712  const int level = -1,
713  const int index = -1,
714  const AccessorData * local_data = nullptr);
715 
728  template <int structdim2, int dim2, int spacedim2>
730 
735  template <int structdim2, int dim2, int spacedim2>
737 
747  void
748  operator=(const TriaAccessor &) = delete;
749 
756  bool
757  used() const;
758 
771  vertex_iterator(const unsigned int i) const;
772 
788  unsigned int
789  vertex_index(const unsigned int i) const;
790 
829  vertex(const unsigned int i) const;
830 
834  typename ::internal::TriangulationImplementation::
835  Iterators<dim, spacedim>::line_iterator
836  line(const unsigned int i) const;
837 
844  unsigned int
845  line_index(const unsigned int i) const;
846 
850  typename ::internal::TriangulationImplementation::
851  Iterators<dim, spacedim>::quad_iterator
852  quad(const unsigned int i) const;
853 
860  unsigned int
861  quad_index(const unsigned int i) const;
884  bool
885  face_orientation(const unsigned int face) const;
886 
896  bool
897  face_flip(const unsigned int face) const;
898 
908  bool
909  face_rotation(const unsigned int face) const;
910 
921  bool
922  line_orientation(const unsigned int line) const;
937  bool
938  has_children() const;
939 
944  unsigned int
945  n_children() const;
946 
960  unsigned int
961  number_of_children() const;
962 
976  unsigned int
977  max_refinement_depth() const;
978 
983  child(const unsigned int i) const;
984 
994  isotropic_child(const unsigned int i) const;
995 
1000  refinement_case() const;
1001 
1007  int
1008  child_index(const unsigned int i) const;
1009 
1015  int
1016  isotropic_child_index(const unsigned int i) const;
1039  boundary_id() const;
1040 
1070  void
1071  set_boundary_id(const types::boundary_id) const;
1072 
1103  void
1104  set_all_boundary_ids(const types::boundary_id) const;
1105 
1113  bool
1114  at_boundary() const;
1115 
1125  const Manifold<dim, spacedim> &
1126  get_manifold() const;
1127 
1149  manifold_id() const;
1150 
1168  void
1169  set_manifold_id(const types::manifold_id) const;
1170 
1184  void
1185  set_all_manifold_ids(const types::manifold_id) const;
1186 
1203  bool
1204  user_flag_set() const;
1205 
1211  void
1212  set_user_flag() const;
1213 
1219  void
1220  clear_user_flag() const;
1221 
1227  void
1228  recursively_set_user_flag() const;
1229 
1235  void
1236  recursively_clear_user_flag() const;
1237 
1243  void
1244  clear_user_data() const;
1245 
1257  void
1258  set_user_pointer(void *p) const;
1259 
1265  void
1266  clear_user_pointer() const;
1267 
1283  void *
1284  user_pointer() const;
1285 
1307  void
1308  recursively_set_user_pointer(void *p) const;
1309 
1316  void
1317  recursively_clear_user_pointer() const;
1318 
1328  void
1329  set_user_index(const unsigned int p) const;
1330 
1336  void
1337  clear_user_index() const;
1338 
1350  unsigned int
1351  user_index() const;
1352 
1370  void
1371  recursively_set_user_index(const unsigned int p) const;
1372 
1381  void
1382  recursively_clear_user_index() const;
1401  double
1402  diameter() const;
1403 
1430  std::pair<Point<spacedim>, double>
1431  enclosing_ball() const;
1432 
1442  bounding_box() const;
1443 
1453  double
1454  extent_in_direction(const unsigned int axis) const;
1455 
1459  double
1460  minimum_vertex_distance() const;
1461 
1476  intermediate_point(const Point<structdim> &coordinates) const;
1477 
1501  real_to_unit_cell_affine_approximation(const Point<spacedim> &point) const;
1502 
1538  center(const bool respect_manifold = false,
1539  const bool interpolate_from_surrounding = false) const;
1540 
1559  barycenter() const;
1560 
1586  double
1587  measure() const;
1588 
1603  bool
1604  is_translation_of(
1606 
1612 private:
1617  void
1618  set_boundary_id_internal(const types::boundary_id id) const;
1619 
1624  void
1625  set(const ::internal::TriangulationImplementation::TriaObject<structdim>
1626  &o) const;
1627 
1635  void
1636  set_line_orientation(const unsigned int line, const bool orientation) const;
1637 
1648  void
1649  set_face_orientation(const unsigned int face, const bool orientation) const;
1650 
1657  void
1658  set_face_flip(const unsigned int face, const bool flip) const;
1659 
1666  void
1667  set_face_rotation(const unsigned int face, const bool rotation) const;
1668 
1672  void
1673  set_used_flag() const;
1674 
1678  void
1679  clear_used_flag() const;
1680 
1689  void
1690  set_refinement_case(const RefinementCase<structdim> &ref_case) const;
1691 
1699  void
1700  clear_refinement_case() const;
1701 
1708  void
1709  set_children(const unsigned int i, const int index) const;
1710 
1715  void
1716  clear_children() const;
1717 
1718 private:
1719  template <int, int>
1720  friend class Triangulation;
1721 
1722  friend struct ::internal::TriangulationImplementation::Implementation;
1723  friend struct ::internal::TriaAccessorImplementation::Implementation;
1724 };
1725 
1726 
1727 
1747 template <int dim, int spacedim>
1748 class TriaAccessor<0, dim, spacedim>
1749 {
1750 public:
1756  static const unsigned int space_dimension = spacedim;
1757 
1763  static const unsigned int dimension = dim;
1764 
1770  static const unsigned int structure_dimension = 0;
1771 
1775  using AccessorData = void;
1776 
1782  const unsigned int vertex_index);
1783 
1789  TriaAccessor(const Triangulation<dim, spacedim> *tria = nullptr,
1790  const int level = 0,
1791  const int index = 0,
1792  const AccessorData * = nullptr);
1793 
1797  template <int structdim2, int dim2, int spacedim2>
1799 
1803  template <int structdim2, int dim2, int spacedim2>
1805 
1810  state() const;
1811 
1816  static int
1817  level();
1818 
1823  int
1824  index() const;
1825 
1831  get_triangulation() const;
1832 
1842  void
1843  operator++();
1844 
1848  void
1849  operator--();
1853  bool
1854  operator==(const TriaAccessor &) const;
1855 
1859  bool
1860  operator!=(const TriaAccessor &) const;
1861 
1889  unsigned int
1890  vertex_index(const unsigned int i = 0) const;
1891 
1897  Point<spacedim> &
1898  vertex(const unsigned int i = 0) const;
1899 
1904  typename ::internal::TriangulationImplementation::
1905  Iterators<dim, spacedim>::line_iterator static line(const unsigned int);
1906 
1910  static unsigned int
1911  line_index(const unsigned int i);
1912 
1916  static typename ::internal::TriangulationImplementation::
1917  Iterators<dim, spacedim>::quad_iterator
1918  quad(const unsigned int i);
1919 
1923  static unsigned int
1924  quad_index(const unsigned int i);
1925 
1941  double
1942  diameter() const;
1943 
1951  double
1952  extent_in_direction(const unsigned int axis) const;
1953 
1962  center(const bool respect_manifold = false,
1963  const bool interpolate_from_surrounding = false) const;
1964 
1973  double
1974  measure() const;
1989  static bool
1990  face_orientation(const unsigned int face);
1991 
1995  static bool
1996  face_flip(const unsigned int face);
1997 
2001  static bool
2002  face_rotation(const unsigned int face);
2003 
2007  static bool
2008  line_orientation(const unsigned int line);
2009 
2024  static bool
2025  has_children();
2026 
2031  static unsigned int
2032  n_children();
2033 
2038  static unsigned int
2039  number_of_children();
2040 
2044  static unsigned int
2045  max_refinement_depth();
2046 
2051  child(const unsigned int);
2052 
2057  isotropic_child(const unsigned int);
2058 
2062  static RefinementCase<0>
2063  refinement_case();
2064 
2068  static int
2069  child_index(const unsigned int i);
2070 
2074  static int
2075  isotropic_child_index(const unsigned int i);
2083  bool
2084  used() const;
2085 
2086 protected:
2094  void
2095  copy_from(const TriaAccessor &);
2096 
2105  bool
2106  operator<(const TriaAccessor &other) const;
2107 
2112 
2116  unsigned int global_vertex_index;
2117 
2118 private:
2119  template <typename Accessor>
2120  friend class TriaRawIterator;
2121  template <typename Accessor>
2122  friend class TriaIterator;
2123  template <typename Accessor>
2124  friend class TriaActiveIterator;
2125 };
2126 
2127 
2128 
2146 template <int spacedim>
2147 class TriaAccessor<0, 1, spacedim>
2148 {
2149 public:
2155  static const unsigned int space_dimension = spacedim;
2156 
2162  static const unsigned int dimension = 1;
2163 
2169  static const unsigned int structure_dimension = 0;
2170 
2174  using AccessorData = void;
2175 
2181  {
2193  right_vertex
2194  };
2195 
2208  const VertexKind vertex_kind,
2209  const unsigned int vertex_index);
2210 
2216  TriaAccessor(const Triangulation<1, spacedim> *tria = nullptr,
2217  const int = 0,
2218  const int = 0,
2219  const AccessorData * = nullptr);
2220 
2224  template <int structdim2, int dim2, int spacedim2>
2226 
2230  template <int structdim2, int dim2, int spacedim2>
2232 
2237  void
2238  copy_from(const TriaAccessor &);
2239 
2246  state();
2247 
2252  static int
2253  level();
2254 
2259  int
2260  index() const;
2261 
2267  get_triangulation() const;
2268 
2279  void
2280  operator++() const;
2281 
2286  void
2287  operator--() const;
2291  bool
2292  operator==(const TriaAccessor &) const;
2293 
2297  bool
2298  operator!=(const TriaAccessor &) const;
2299 
2308  bool
2309  operator<(const TriaAccessor &other) const;
2310 
2337  unsigned int
2338  vertex_index(const unsigned int i = 0) const;
2339 
2345  Point<spacedim> &
2346  vertex(const unsigned int i = 0) const;
2347 
2353  center() const;
2354 
2359  typename ::internal::TriangulationImplementation::
2360  Iterators<1, spacedim>::line_iterator static line(const unsigned int);
2361 
2368  static unsigned int
2369  line_index(const unsigned int i);
2370 
2374  static typename ::internal::TriangulationImplementation::
2375  Iterators<1, spacedim>::quad_iterator
2376  quad(const unsigned int i);
2377 
2384  static unsigned int
2385  quad_index(const unsigned int i);
2386 
2396  bool
2397  at_boundary() const;
2398 
2414  boundary_id() const;
2415 
2419  const Manifold<1, spacedim> &
2420  get_manifold() const;
2421 
2429  manifold_id() const;
2430 
2431 
2442  static bool
2443  face_orientation(const unsigned int face);
2444 
2448  static bool
2449  face_flip(const unsigned int face);
2450 
2454  static bool
2455  face_rotation(const unsigned int face);
2456 
2460  static bool
2461  line_orientation(const unsigned int line);
2462 
2477  static bool
2478  has_children();
2479 
2484  static unsigned int
2485  n_children();
2486 
2491  static unsigned int
2492  number_of_children();
2493 
2497  static unsigned int
2498  max_refinement_depth();
2499 
2504  child(const unsigned int);
2505 
2510  isotropic_child(const unsigned int);
2511 
2515  static RefinementCase<0>
2516  refinement_case();
2517 
2521  static int
2522  child_index(const unsigned int i);
2523 
2527  static int
2528  isotropic_child_index(const unsigned int i);
2561  void
2562  set_boundary_id(const types::boundary_id);
2563 
2570  void
2571  set_manifold_id(const types::manifold_id);
2572 
2584  void
2585  set_all_boundary_ids(const types::boundary_id);
2586 
2598  void
2599  set_all_manifold_ids(const types::manifold_id);
2607  bool
2608  used() const;
2609 
2610 protected:
2615 
2621 
2625  unsigned int global_vertex_index;
2626 };
2627 
2628 
2629 
2646 template <int dim, int spacedim = dim>
2647 class CellAccessor : public TriaAccessor<dim, dim, spacedim>
2648 {
2649 public:
2654 
2659 
2670  CellAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
2671  const int level = -1,
2672  const int index = -1,
2673  const AccessorData * local_data = nullptr);
2674 
2678  CellAccessor(const TriaAccessor<dim, dim, spacedim> &cell_accessor);
2679 
2692  template <int structdim2, int dim2, int spacedim2>
2694 
2699  template <int structdim2, int dim2, int spacedim2>
2701 
2711  void
2712  operator=(const CellAccessor<dim, spacedim> &) = delete;
2713 
2730  child(const unsigned int i) const;
2731 
2735  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>
2736  face(const unsigned int i) const;
2737 
2738 
2742  std::array<TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>,
2744  face_iterators() const;
2745 
2755  unsigned int
2756  face_index(const unsigned int i) const;
2757 
2806  neighbor_child_on_subface(const unsigned int face_no,
2807  const unsigned int subface_no) const;
2808 
2831  neighbor(const unsigned int i) const;
2832 
2837  int
2838  neighbor_index(const unsigned int i) const;
2839 
2844  int
2845  neighbor_level(const unsigned int i) const;
2846 
2858  unsigned int
2859  neighbor_of_neighbor(const unsigned int neighbor) const;
2860 
2871  bool
2872  neighbor_is_coarser(const unsigned int neighbor) const;
2873 
2888  std::pair<unsigned int, unsigned int>
2889  neighbor_of_coarser_neighbor(const unsigned int neighbor) const;
2890 
2897  unsigned int
2898  neighbor_face_no(const unsigned int neighbor) const;
2899 
2903  static bool
2904  is_level_cell();
2905 
2919  bool
2920  has_periodic_neighbor(const unsigned int i) const;
2921 
2939  periodic_neighbor(const unsigned int i) const;
2940 
2949  neighbor_or_periodic_neighbor(const unsigned int i) const;
2950 
2966  periodic_neighbor_child_on_subface(const unsigned int face_no,
2967  const unsigned int subface_no) const;
2968 
2979  std::pair<unsigned int, unsigned int>
2980  periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const;
2981 
2987  int
2988  periodic_neighbor_index(const unsigned int i) const;
2989 
2995  int
2996  periodic_neighbor_level(const unsigned int i) const;
2997 
3012  unsigned int
3013  periodic_neighbor_of_periodic_neighbor(const unsigned int i) const;
3014 
3020  unsigned int
3021  periodic_neighbor_face_no(const unsigned int i) const;
3022 
3029  bool
3030  periodic_neighbor_is_coarser(const unsigned int i) const;
3031 
3048  bool
3049  at_boundary(const unsigned int i) const;
3050 
3059  bool
3060  at_boundary() const;
3061 
3069  bool
3070  has_boundary_lines() const;
3097  refine_flag_set() const;
3098 
3116  void
3117  set_refine_flag(const RefinementCase<dim> ref_case =
3119 
3123  void
3124  clear_refine_flag() const;
3125 
3133  bool
3134  flag_for_face_refinement(
3135  const unsigned int face_no,
3136  const RefinementCase<dim - 1> &face_refinement_case =
3138 
3144  bool
3145  flag_for_line_refinement(const unsigned int line_no) const;
3146 
3156  subface_case(const unsigned int face_no) const;
3157 
3161  bool
3162  coarsen_flag_set() const;
3163 
3168  void
3169  set_coarsen_flag() const;
3170 
3174  void
3175  clear_coarsen_flag() const;
3199  material_id() const;
3200 
3212  void
3213  set_material_id(const types::material_id new_material_id) const;
3214 
3223  void
3224  recursively_set_material_id(const types::material_id new_material_id) const;
3251  subdomain_id() const;
3252 
3268  void
3269  set_subdomain_id(const types::subdomain_id new_subdomain_id) const;
3270 
3280  level_subdomain_id() const;
3281 
3286  void
3287  set_level_subdomain_id(
3288  const types::subdomain_id new_level_subdomain_id) const;
3289 
3290 
3306  void
3307  recursively_set_subdomain_id(
3308  const types::subdomain_id new_subdomain_id) const;
3326  bool
3327  direction_flag() const;
3328 
3354  unsigned int
3355  active_cell_index() const;
3356 
3364  int
3365  parent_index() const;
3366 
3373  parent() const;
3374 
3394  bool
3395  active() const;
3396 
3416  bool
3417  is_locally_owned() const;
3418 
3423  bool
3424  is_locally_owned_on_level() const;
3425 
3449  bool
3450  is_ghost() const;
3451 
3478  bool
3479  is_artificial() const;
3480 
3494  bool
3495  point_inside(const Point<spacedim> &p) const;
3496 
3505  void
3506  set_neighbor(const unsigned int i,
3507  const TriaIterator<CellAccessor<dim, spacedim>> &pointer) const;
3508 
3522  CellId
3523  id() const;
3524 
3533  DeclException0(ExcRefineCellNotActive);
3537  DeclException0(ExcCellFlaggedForRefinement);
3541  DeclException0(ExcCellFlaggedForCoarsening);
3542 
3543 protected:
3559  unsigned int
3560  neighbor_of_neighbor_internal(const unsigned int neighbor) const;
3561 
3567  template <int dim_, int spacedim_>
3568  bool
3569  point_inside_codim(const Point<spacedim_> &p) const;
3570 
3571 
3572 
3573 private:
3578  void
3579  set_active_cell_index(const unsigned int active_cell_index);
3580 
3584  void
3585  set_parent(const unsigned int parent_index);
3586 
3593  void
3594  set_direction_flag(const bool new_direction_flag) const;
3595 
3596  template <int, int>
3597  friend class Triangulation;
3598 
3599  friend struct ::internal::TriangulationImplementation::Implementation;
3600 };
3601 
3602 
3603 
3604 /* ----- declaration of explicit specializations and general templates ----- */
3605 
3606 
3607 template <int structdim, int dim, int spacedim>
3608 template <typename OtherAccessor>
3610  const OtherAccessor &)
3611 {
3612  Assert(false,
3613  ExcMessage("You are attempting an illegal conversion between "
3614  "iterator/accessor types. The constructor you call "
3615  "only exists to make certain template constructs "
3616  "easier to write as dimension independent code but "
3617  "the conversion is not valid in the current context."));
3618 }
3619 
3620 
3621 
3622 template <int structdim, int dim, int spacedim>
3623 template <int structdim2, int dim2, int spacedim2>
3626 {
3627  Assert(false,
3628  ExcMessage("You are attempting an illegal conversion between "
3629  "iterator/accessor types. The constructor you call "
3630  "only exists to make certain template constructs "
3631  "easier to write as dimension independent code but "
3632  "the conversion is not valid in the current context."));
3633 }
3634 
3635 
3636 
3637 template <int dim, int spacedim>
3638 template <int structdim2, int dim2, int spacedim2>
3641 {
3642  Assert(false,
3643  ExcMessage("You are attempting an illegal conversion between "
3644  "iterator/accessor types. The constructor you call "
3645  "only exists to make certain template constructs "
3646  "easier to write as dimension independent code but "
3647  "the conversion is not valid in the current context."));
3648 }
3649 
3650 
3651 
3652 template <int structdim, int dim, int spacedim>
3653 template <int structdim2, int dim2, int spacedim2>
3656 {
3657  Assert(false,
3658  ExcMessage("You are attempting an illegal conversion between "
3659  "iterator/accessor types. The constructor you call "
3660  "only exists to make certain template constructs "
3661  "easier to write as dimension independent code but "
3662  "the conversion is not valid in the current context."));
3663 }
3664 
3665 
3666 
3667 template <int dim, int spacedim>
3668 template <int structdim2, int dim2, int spacedim2>
3671 {
3672  Assert(false,
3673  ExcMessage("You are attempting an illegal conversion between "
3674  "iterator/accessor types. The constructor you call "
3675  "only exists to make certain template constructs "
3676  "easier to write as dimension independent code but "
3677  "the conversion is not valid in the current context."));
3678 }
3679 
3680 
3681 #ifndef DOXYGEN
3682 
3683 template <>
3684 bool
3686 template <>
3687 bool
3689 template <>
3690 bool
3692 template <>
3693 bool
3695 template <>
3696 bool
3698 template <>
3699 bool
3701 // -------------------------------------------------------------------
3702 
3703 template <>
3704 void
3706 
3707 #endif // DOXYGEN
3708 
3709 DEAL_II_NAMESPACE_CLOSE
3710 
3711 // include more templates in debug and optimized mode
3712 #include "tria_accessor.templates.h"
3713 
3714 #endif
unsigned int manifold_id
Definition: types.h:137
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()
const Triangulation< dim, spacedim > * tria
static ::ExceptionBase & ExcNeighborIsCoarser()
static ::ExceptionBase & ExcDereferenceInvalidObject(AccessorType arg1)
const Triangulation< 1, spacedim > * tria
unsigned int material_id
Definition: types.h:148
void set_all_manifold_ids(const types::manifold_id) const
static ::ExceptionBase & ExcFacesHaveNoLevel()
static ::ExceptionBase & ExcCellNotUsed()
static ::ExceptionBase & ExcMessage(std::string arg1)
unsigned int subdomain_id
Definition: types.h:43
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:518
const Triangulation< dim, spacedim > * tria
static ::ExceptionBase & ExcSetOnlyEvenChildren(int arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1407
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:496
#define DeclException0(Exception0)
Definition: exceptions.h:473
static ::ExceptionBase & ExcCellNotActive()
static ::ExceptionBase & ExcNoPeriodicNeighbor()
static ::ExceptionBase & ExcCantSetChildren(int arg1)
Definition: cell_id.h:68
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
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.
unsigned int boundary_id
Definition: types.h:125
typename ::internal::TriaAccessorImplementation::PresentLevelType< structdim, dim >::type present_level
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcCantCompareIterators()