Reference documentation for deal.II version Git 409ee4b167 2020-08-14 09:46:12 -0400
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
tria_accessor.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2020 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_tria_accessor_h
17 #define dealii_tria_accessor_h
18 
19 
20 #include <deal.II/base/config.h>
21 
25 #include <deal.II/base/point.h>
26 
27 #include <deal.II/grid/cell_id.h>
31 
32 #include <utility>
33 
34 
36 
37 // Forward declarations
38 #ifndef DOXYGEN
39 template <int dim, int spacedim>
40 class Triangulation;
41 template <typename Accessor>
42 class TriaRawIterator;
43 template <typename Accessor>
44 class TriaIterator;
45 template <typename Accessor>
46 class TriaActiveIterator;
47 
48 template <int dim, int spacedim>
49 class Manifold;
50 #endif
51 
52 namespace internal
53 {
54  namespace TriangulationImplementation
55  {
56  class TriaObjects;
57  struct Implementation;
58  } // namespace TriangulationImplementation
59 
60  namespace TriaAccessorImplementation
61  {
62  struct Implementation;
63 
69  template <int structdim, int dim>
71  {
72  struct type
73  {
77  type() = default;
78 
82  type(const int level)
83  {
84  Assert(level == 0, ExcInternalError());
85  (void)level; // removes -Wunused-parameter warning in optimized mode
86  }
87 
91  operator int() const
92  {
93  return 0;
94  }
95 
96  void
97  operator++() const
98  {
99  Assert(false, ExcInternalError());
100  }
101 
102  void
103  operator--() const
104  {
105  Assert(false, ExcInternalError());
106  }
107  };
108  };
109 
110 
116  template <int dim>
117  struct PresentLevelType<dim, dim>
118  {
119  using type = int;
120  };
121 
122  } // namespace TriaAccessorImplementation
123 } // namespace internal
124 template <int structdim, int dim, int spacedim>
126 template <int dim, int spacedim>
127 class TriaAccessor<0, dim, spacedim>;
128 template <int spacedim>
129 class TriaAccessor<0, 1, spacedim>;
130 
131 // note: the file tria_accessor.templates.h is included at the end of
132 // this file. this includes a lot of templates. originally, this was
133 // only done in debug mode, but led to cyclic reduction problems and
134 // so is now on by default.
135 
136 
141 {
146  "The operation you are attempting can only be performed for "
147  "(cell, face, or edge) iterators that point to valid "
148  "objects. These objects need not necessarily be active, "
149  "i.e., have no children, but they need to be part of a "
150  "triangulation. (The objects pointed to by an iterator "
151  "may -- after coarsening -- also be objects that used "
152  "to be part of a triangulation, but are now no longer "
153  "used. Their memory location may have been retained "
154  "for re-use upon the next mesh refinement, but is "
155  "currently unused.)");
166  "The operation you are attempting can only be performed for "
167  "(cell, face, or edge) iterators that point to 'active' "
168  "objects. 'Active' objects are those that do not have "
169  "children (in the case of cells), or that are part of "
170  "an active cell (in the case of faces or edges). However, "
171  "the object on which you are trying the current "
172  "operation is not 'active' in this sense.");
179  "The operation you are attempting can only be performed for "
180  "(cell, face, or edge) iterators that have children, "
181  "but the object on which you are trying the current "
182  "operation does not have any.");
190  "The operation you are attempting can only be performed for "
191  "(cell, face, or edge) iterators that have a parent object, "
192  "but the object on which you are trying the current "
193  "operation does not have one -- i.e., it is on the "
194  "coarsest level of the triangulation.");
199  int,
200  << "You can only set the child index if the cell does not "
201  << "currently have children registered; or you can clear it. "
202  << "The given index was " << arg1
203  << " (-1 means: clear children).");
207  template <typename AccessorType>
209  AccessorType,
210  << "You tried to dereference an iterator for which this "
211  << "is not possible. More information on this iterator: "
212  << "index=" << arg1.index() << ", state="
213  << (arg1.state() == IteratorState::valid ?
214  "valid" :
215  (arg1.state() == IteratorState::past_the_end ?
216  "past_the_end" :
217  "invalid")));
222  "Iterators can only be compared if they point to the same "
223  "triangulation, or if neither of them are associated "
224  "with a triangulation.");
225  // TODO: Write documentation!
230  // TODO: Write documentation!
250  // TODO: Write documentation!
256  int,
257  << "You can only set the child index of an even numbered child."
258  << "The number of the child given was " << arg1 << ".");
259 } // namespace TriaAccessorExceptions
260 
261 
287 template <int structdim, int dim, int spacedim = dim>
289 {
290 public:
296  static const unsigned int space_dimension = spacedim;
297 
303  static const unsigned int dimension = dim;
304 
310  static const unsigned int structure_dimension = structdim;
311 
321  void
322  operator=(const TriaAccessorBase *) = delete;
323 
324 protected:
330  using AccessorData = void;
331 
335  TriaAccessorBase(const Triangulation<dim, spacedim> *parent = nullptr,
336  const int level = -1,
337  const int index = -1,
338  const AccessorData * = nullptr);
339 
344 
352  void
353  copy_from(const TriaAccessorBase &);
354 
359  operator=(const TriaAccessorBase &);
360 
371  bool
372  operator<(const TriaAccessorBase &other) const;
373 
374 protected:
378  bool
379  operator==(const TriaAccessorBase &) const;
380 
384  bool
385  operator!=(const TriaAccessorBase &) const;
386 
400  void
401  operator++();
402 
410  void
411  operator--();
420  objects() const;
421 
422 public:
428  using LocalData = void *;
429 
453  int
454  level() const;
455 
482  int
483  index() const;
484 
490  state() const;
491 
497  get_triangulation() const;
498 
502 protected:
507  typename ::internal::TriaAccessorImplementation::
508  PresentLevelType<structdim, dim>::type present_level;
509 
515 
520 
521 private:
522  template <typename Accessor>
523  friend class TriaRawIterator;
524  template <typename Accessor>
525  friend class TriaIterator;
526  template <typename Accessor>
527  friend class TriaActiveIterator;
528 };
529 
530 
531 
552 template <int structdim, int dim, int spacedim = dim>
553 class InvalidAccessor : public TriaAccessorBase<structdim, dim, spacedim>
554 {
555 public:
559  using AccessorData =
561 
569  InvalidAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
570  const int level = -1,
571  const int index = -1,
572  const AccessorData * local_data = nullptr);
573 
582 
587  template <typename OtherAccessor>
588  InvalidAccessor(const OtherAccessor &);
589 
593  void
594  copy_from(const InvalidAccessor &);
595 
599  bool
600  operator==(const InvalidAccessor &) const;
601  bool
602  operator!=(const InvalidAccessor &) const;
603 
607  void
608  operator++() const;
609  void
610  operator--() const;
611 
616  bool
617  used() const;
618 
623  bool
624  has_children() const;
625 
630  manifold_id() const;
631 
635  unsigned int
636  user_index() const;
637 
641  void
642  set_user_index(const unsigned int p) const;
643 
647  void
648  set_manifold_id(const types::manifold_id) const;
649 
654  vertex(const unsigned int i) const;
655 
660  typename ::internal::TriangulationImplementation::
661  Iterators<dim, spacedim>::line_iterator
662  line(const unsigned int i) const;
663 
668  typename ::internal::TriangulationImplementation::
669  Iterators<dim, spacedim>::quad_iterator
670  quad(const unsigned int i) const;
671 };
672 
673 
674 
692 template <int structdim, int dim, int spacedim>
693 class TriaAccessor : public TriaAccessorBase<structdim, dim, spacedim>
694 {
695 public:
699  using AccessorData =
701 
705  TriaAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
706  const int level = -1,
707  const int index = -1,
708  const AccessorData * local_data = nullptr);
709 
722  template <int structdim2, int dim2, int spacedim2>
724 
729  template <int structdim2, int dim2, int spacedim2>
731 
741  void
742  operator=(const TriaAccessor &) = delete;
743 
750  bool
751  used() const;
752 
765  vertex_iterator(const unsigned int i) const;
766 
782  unsigned int
783  vertex_index(const unsigned int i) const;
784 
823  vertex(const unsigned int i) const;
824 
828  typename ::internal::TriangulationImplementation::
829  Iterators<dim, spacedim>::line_iterator
830  line(const unsigned int i) const;
831 
838  unsigned int
839  line_index(const unsigned int i) const;
840 
844  typename ::internal::TriangulationImplementation::
845  Iterators<dim, spacedim>::quad_iterator
846  quad(const unsigned int i) const;
847 
854  unsigned int
855  quad_index(const unsigned int i) const;
878  bool
879  face_orientation(const unsigned int face) const;
880 
890  bool
891  face_flip(const unsigned int face) const;
892 
902  bool
903  face_rotation(const unsigned int face) const;
904 
915  bool
916  line_orientation(const unsigned int line) const;
931  bool
932  has_children() const;
933 
938  unsigned int
939  n_children() const;
940 
954  unsigned int
955  number_of_children() const;
956 
970  unsigned int
971  max_refinement_depth() const;
972 
977  child(const unsigned int i) const;
978 
983  unsigned int
984  child_iterator_to_index(
986 
996  isotropic_child(const unsigned int i) const;
997 
1002  refinement_case() const;
1003 
1009  int
1010  child_index(const unsigned int i) const;
1011 
1017  int
1018  isotropic_child_index(const unsigned int i) const;
1041  boundary_id() const;
1042 
1072  void
1073  set_boundary_id(const types::boundary_id) const;
1074 
1105  void
1106  set_all_boundary_ids(const types::boundary_id) const;
1107 
1115  bool
1116  at_boundary() const;
1117 
1127  const Manifold<dim, spacedim> &
1128  get_manifold() const;
1129 
1151  manifold_id() const;
1152 
1170  void
1171  set_manifold_id(const types::manifold_id) const;
1172 
1186  void
1187  set_all_manifold_ids(const types::manifold_id) const;
1188 
1205  bool
1206  user_flag_set() const;
1207 
1213  void
1214  set_user_flag() const;
1215 
1221  void
1222  clear_user_flag() const;
1223 
1229  void
1230  recursively_set_user_flag() const;
1231 
1237  void
1238  recursively_clear_user_flag() const;
1239 
1245  void
1246  clear_user_data() const;
1247 
1259  void
1260  set_user_pointer(void *p) const;
1261 
1267  void
1268  clear_user_pointer() const;
1269 
1285  void *
1286  user_pointer() const;
1287 
1309  void
1310  recursively_set_user_pointer(void *p) const;
1311 
1318  void
1319  recursively_clear_user_pointer() const;
1320 
1330  void
1331  set_user_index(const unsigned int p) const;
1332 
1338  void
1339  clear_user_index() const;
1340 
1352  unsigned int
1353  user_index() const;
1354 
1372  void
1373  recursively_set_user_index(const unsigned int p) const;
1374 
1383  void
1384  recursively_clear_user_index() const;
1403  double
1404  diameter() const;
1405 
1432  std::pair<Point<spacedim>, double>
1433  enclosing_ball() const;
1434 
1444  bounding_box() const;
1445 
1455  double
1456  extent_in_direction(const unsigned int axis) const;
1457 
1461  double
1462  minimum_vertex_distance() const;
1463 
1478  intermediate_point(const Point<structdim> &coordinates) const;
1479 
1503  real_to_unit_cell_affine_approximation(const Point<spacedim> &point) const;
1504 
1540  center(const bool respect_manifold = false,
1541  const bool interpolate_from_surrounding = false) const;
1542 
1561  barycenter() const;
1562 
1588  double
1589  measure() const;
1590 
1605  bool
1606  is_translation_of(
1608 
1613  reference_cell_type() const;
1614 
1618  unsigned int
1619  n_vertices() const;
1620 
1624  unsigned int
1625  n_lines() const;
1626 
1632  unsigned int
1633  n_faces() const;
1634 
1640  vertex_indices() const;
1641 
1647  line_indices() const;
1648 
1656  face_indices() const;
1657 
1662 protected:
1668  reference_cell_info() const;
1669 
1670 private:
1675  void
1676  set_boundary_id_internal(const types::boundary_id id) const;
1677 
1685  void
1686  set_bounding_object_indices(
1687  const std::initializer_list<int> &new_indices) const;
1688 
1692  void
1693  set_bounding_object_indices(
1694  const std::initializer_list<unsigned int> &new_indices) const;
1695 
1703  void
1704  set_line_orientation(const unsigned int line, const bool orientation) const;
1705 
1716  void
1717  set_face_orientation(const unsigned int face, const bool orientation) const;
1718 
1725  void
1726  set_face_flip(const unsigned int face, const bool flip) const;
1727 
1734  void
1735  set_face_rotation(const unsigned int face, const bool rotation) const;
1736 
1740  void
1741  set_used_flag() const;
1742 
1746  void
1747  clear_used_flag() const;
1748 
1757  void
1758  set_refinement_case(const RefinementCase<structdim> &ref_case) const;
1759 
1767  void
1768  clear_refinement_case() const;
1769 
1776  void
1777  set_children(const unsigned int i, const int index) const;
1778 
1783  void
1784  clear_children() const;
1785 
1786 private:
1787  template <int, int>
1788  friend class Triangulation;
1789 
1790  friend struct ::internal::TriangulationImplementation::Implementation;
1791  friend struct ::internal::TriaAccessorImplementation::Implementation;
1792 };
1793 
1794 
1795 
1814 template <int dim, int spacedim>
1815 class TriaAccessor<0, dim, spacedim>
1816 {
1817 public:
1823  static const unsigned int space_dimension = spacedim;
1824 
1830  static const unsigned int dimension = dim;
1831 
1837  static const unsigned int structure_dimension = 0;
1838 
1842  using AccessorData = void;
1843 
1849  const unsigned int vertex_index);
1850 
1856  TriaAccessor(const Triangulation<dim, spacedim> *tria = nullptr,
1857  const int level = 0,
1858  const int index = 0,
1859  const AccessorData * = nullptr);
1860 
1864  template <int structdim2, int dim2, int spacedim2>
1866 
1870  template <int structdim2, int dim2, int spacedim2>
1872 
1877  state() const;
1878 
1883  static int
1884  level();
1885 
1890  int
1891  index() const;
1892 
1898  get_triangulation() const;
1899 
1909  void
1910  operator++();
1911 
1915  void
1916  operator--();
1920  bool
1921  operator==(const TriaAccessor &) const;
1922 
1926  bool
1927  operator!=(const TriaAccessor &) const;
1928 
1956  unsigned int
1957  vertex_index(const unsigned int i = 0) const;
1958 
1964  Point<spacedim> &
1965  vertex(const unsigned int i = 0) const;
1966 
1971  typename ::internal::TriangulationImplementation::
1972  Iterators<dim, spacedim>::line_iterator static line(const unsigned int);
1973 
1977  static unsigned int
1978  line_index(const unsigned int i);
1979 
1983  static typename ::internal::TriangulationImplementation::
1984  Iterators<dim, spacedim>::quad_iterator
1985  quad(const unsigned int i);
1986 
1990  static unsigned int
1991  quad_index(const unsigned int i);
1992 
2008  double
2009  diameter() const;
2010 
2018  double
2019  extent_in_direction(const unsigned int axis) const;
2020 
2029  center(const bool respect_manifold = false,
2030  const bool interpolate_from_surrounding = false) const;
2031 
2040  double
2041  measure() const;
2056  static bool
2057  face_orientation(const unsigned int face);
2058 
2062  static bool
2063  face_flip(const unsigned int face);
2064 
2068  static bool
2069  face_rotation(const unsigned int face);
2070 
2074  static bool
2075  line_orientation(const unsigned int line);
2076 
2091  static bool
2092  has_children();
2093 
2098  static unsigned int
2099  n_children();
2100 
2105  static unsigned int
2106  number_of_children();
2107 
2111  static unsigned int
2112  max_refinement_depth();
2113 
2117  static unsigned int
2118  child_iterator_to_index(const TriaIterator<TriaAccessor<0, dim, spacedim>> &);
2119 
2124  child(const unsigned int);
2125 
2130  isotropic_child(const unsigned int);
2131 
2135  static RefinementCase<0>
2136  refinement_case();
2137 
2141  static int
2142  child_index(const unsigned int i);
2143 
2147  static int
2148  isotropic_child_index(const unsigned int i);
2156  bool
2157  used() const;
2158 
2159 protected:
2167  void
2168  copy_from(const TriaAccessor &);
2169 
2178  bool
2179  operator<(const TriaAccessor &other) const;
2180 
2185 
2189  unsigned int global_vertex_index;
2190 
2191 private:
2192  template <typename Accessor>
2193  friend class TriaRawIterator;
2194  template <typename Accessor>
2195  friend class TriaIterator;
2196  template <typename Accessor>
2197  friend class TriaActiveIterator;
2198 };
2199 
2200 
2201 
2218 template <int spacedim>
2219 class TriaAccessor<0, 1, spacedim>
2220 {
2221 public:
2227  static const unsigned int space_dimension = spacedim;
2228 
2234  static const unsigned int dimension = 1;
2235 
2241  static const unsigned int structure_dimension = 0;
2242 
2246  using AccessorData = void;
2247 
2253  {
2265  right_vertex
2266  };
2267 
2280  const VertexKind vertex_kind,
2281  const unsigned int vertex_index);
2282 
2288  TriaAccessor(const Triangulation<1, spacedim> *tria = nullptr,
2289  const int = 0,
2290  const int = 0,
2291  const AccessorData * = nullptr);
2292 
2296  template <int structdim2, int dim2, int spacedim2>
2298 
2302  template <int structdim2, int dim2, int spacedim2>
2304 
2309  void
2310  copy_from(const TriaAccessor &);
2311 
2318  state();
2319 
2324  static int
2325  level();
2326 
2331  int
2332  index() const;
2333 
2339  get_triangulation() const;
2340 
2351  void
2352  operator++() const;
2353 
2358  void
2359  operator--() const;
2363  bool
2364  operator==(const TriaAccessor &) const;
2365 
2369  bool
2370  operator!=(const TriaAccessor &) const;
2371 
2380  bool
2381  operator<(const TriaAccessor &other) const;
2382 
2409  unsigned int
2410  vertex_index(const unsigned int i = 0) const;
2411 
2417  Point<spacedim> &
2418  vertex(const unsigned int i = 0) const;
2419 
2425  center() const;
2426 
2431  typename ::internal::TriangulationImplementation::
2432  Iterators<1, spacedim>::line_iterator static line(const unsigned int);
2433 
2440  static unsigned int
2441  line_index(const unsigned int i);
2442 
2446  static typename ::internal::TriangulationImplementation::
2447  Iterators<1, spacedim>::quad_iterator
2448  quad(const unsigned int i);
2449 
2456  static unsigned int
2457  quad_index(const unsigned int i);
2458 
2468  bool
2469  at_boundary() const;
2470 
2486  boundary_id() const;
2487 
2491  const Manifold<1, spacedim> &
2492  get_manifold() const;
2493 
2501  manifold_id() const;
2502 
2503 
2514  static bool
2515  face_orientation(const unsigned int face);
2516 
2520  static bool
2521  face_flip(const unsigned int face);
2522 
2526  static bool
2527  face_rotation(const unsigned int face);
2528 
2532  static bool
2533  line_orientation(const unsigned int line);
2534 
2549  static bool
2550  has_children();
2551 
2556  static unsigned int
2557  n_children();
2558 
2563  static unsigned int
2564  number_of_children();
2565 
2569  static unsigned int
2570  max_refinement_depth();
2571 
2575  static unsigned int
2576  child_iterator_to_index(const TriaIterator<TriaAccessor<0, 1, spacedim>> &);
2577 
2582  child(const unsigned int);
2583 
2588  isotropic_child(const unsigned int);
2589 
2593  static RefinementCase<0>
2594  refinement_case();
2595 
2599  static int
2600  child_index(const unsigned int i);
2601 
2605  static int
2606  isotropic_child_index(const unsigned int i);
2639  void
2640  set_boundary_id(const types::boundary_id);
2641 
2648  void
2649  set_manifold_id(const types::manifold_id);
2650 
2662  void
2663  set_all_boundary_ids(const types::boundary_id);
2664 
2676  void
2677  set_all_manifold_ids(const types::manifold_id);
2685  bool
2686  used() const;
2687 
2692  reference_cell_type() const;
2693 
2697  unsigned int
2698  n_vertices() const;
2699 
2703  unsigned int
2704  n_lines() const;
2705 
2711  vertex_indices() const;
2712 
2718  line_indices() const;
2719 
2720 protected:
2725 
2731 
2735  unsigned int global_vertex_index;
2736 };
2737 
2738 
2739 
2755 template <int dim, int spacedim = dim>
2756 class CellAccessor : public TriaAccessor<dim, dim, spacedim>
2757 {
2758 public:
2763 
2768 
2779  CellAccessor(const Triangulation<dim, spacedim> *parent = nullptr,
2780  const int level = -1,
2781  const int index = -1,
2782  const AccessorData * local_data = nullptr);
2783 
2787  CellAccessor(const TriaAccessor<dim, dim, spacedim> &cell_accessor);
2788 
2801  template <int structdim2, int dim2, int spacedim2>
2803 
2808  template <int structdim2, int dim2, int spacedim2>
2810 
2820  void
2821  operator=(const CellAccessor<dim, spacedim> &) = delete;
2822 
2839  child(const unsigned int i) const;
2840 
2844  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>
2845  face(const unsigned int i) const;
2846 
2851  unsigned int
2852  face_iterator_to_index(
2854 
2858  boost::container::small_vector<
2859  TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>,
2861  face_iterators() const;
2862 
2872  unsigned int
2873  face_index(const unsigned int i) const;
2874 
2923  neighbor_child_on_subface(const unsigned int face_no,
2924  const unsigned int subface_no) const;
2925 
2948  neighbor(const unsigned int i) const;
2949 
2954  int
2955  neighbor_index(const unsigned int i) const;
2956 
2961  int
2962  neighbor_level(const unsigned int i) const;
2963 
2975  unsigned int
2976  neighbor_of_neighbor(const unsigned int neighbor) const;
2977 
2988  bool
2989  neighbor_is_coarser(const unsigned int neighbor) const;
2990 
3005  std::pair<unsigned int, unsigned int>
3006  neighbor_of_coarser_neighbor(const unsigned int neighbor) const;
3007 
3014  unsigned int
3015  neighbor_face_no(const unsigned int neighbor) const;
3016 
3020  static bool
3021  is_level_cell();
3022 
3036  bool
3037  has_periodic_neighbor(const unsigned int i) const;
3038 
3056  periodic_neighbor(const unsigned int i) const;
3057 
3066  neighbor_or_periodic_neighbor(const unsigned int i) const;
3067 
3083  periodic_neighbor_child_on_subface(const unsigned int face_no,
3084  const unsigned int subface_no) const;
3085 
3096  std::pair<unsigned int, unsigned int>
3097  periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const;
3098 
3104  int
3105  periodic_neighbor_index(const unsigned int i) const;
3106 
3112  int
3113  periodic_neighbor_level(const unsigned int i) const;
3114 
3129  unsigned int
3130  periodic_neighbor_of_periodic_neighbor(const unsigned int i) const;
3131 
3137  unsigned int
3138  periodic_neighbor_face_no(const unsigned int i) const;
3139 
3146  bool
3147  periodic_neighbor_is_coarser(const unsigned int i) const;
3148 
3165  bool
3166  at_boundary(const unsigned int i) const;
3167 
3176  bool
3177  at_boundary() const;
3178 
3186  bool
3187  has_boundary_lines() const;
3214  refine_flag_set() const;
3215 
3233  void
3234  set_refine_flag(const RefinementCase<dim> ref_case =
3236 
3240  void
3241  clear_refine_flag() const;
3242 
3250  bool
3251  flag_for_face_refinement(
3252  const unsigned int face_no,
3253  const RefinementCase<dim - 1> &face_refinement_case =
3255 
3261  bool
3262  flag_for_line_refinement(const unsigned int line_no) const;
3263 
3273  subface_case(const unsigned int face_no) const;
3274 
3278  bool
3279  coarsen_flag_set() const;
3280 
3285  void
3286  set_coarsen_flag() const;
3287 
3291  void
3292  clear_coarsen_flag() const;
3316  material_id() const;
3317 
3329  void
3330  set_material_id(const types::material_id new_material_id) const;
3331 
3340  void
3341  recursively_set_material_id(const types::material_id new_material_id) const;
3368  subdomain_id() const;
3369 
3385  void
3386  set_subdomain_id(const types::subdomain_id new_subdomain_id) const;
3387 
3397  level_subdomain_id() const;
3398 
3403  void
3404  set_level_subdomain_id(
3405  const types::subdomain_id new_level_subdomain_id) const;
3406 
3407 
3423  void
3424  recursively_set_subdomain_id(
3425  const types::subdomain_id new_subdomain_id) const;
3443  bool
3444  direction_flag() const;
3445 
3471  unsigned int
3472  active_cell_index() const;
3473 
3481  int
3482  parent_index() const;
3483 
3490  parent() const;
3491 
3517  bool
3518  active() const;
3519 
3528  bool
3529  is_active() const;
3530 
3550  bool
3551  is_locally_owned() const;
3552 
3557  bool
3558  is_locally_owned_on_level() const;
3559 
3583  bool
3584  is_ghost() const;
3585 
3612  bool
3613  is_artificial() const;
3614 
3628  bool
3629  point_inside(const Point<spacedim> &p) const;
3630 
3639  void
3640  set_neighbor(const unsigned int i,
3641  const TriaIterator<CellAccessor<dim, spacedim>> &pointer) const;
3642 
3656  CellId
3657  id() const;
3658 
3667  DeclException0(ExcRefineCellNotActive);
3671  DeclException0(ExcCellFlaggedForRefinement);
3675  DeclException0(ExcCellFlaggedForCoarsening);
3676 
3677 protected:
3693  unsigned int
3694  neighbor_of_neighbor_internal(const unsigned int neighbor) const;
3695 
3701  template <int dim_, int spacedim_>
3702  bool
3703  point_inside_codim(const Point<spacedim_> &p) const;
3704 
3705 
3706 
3707 private:
3712  void
3713  set_active_cell_index(const unsigned int active_cell_index);
3714 
3718  void
3719  set_parent(const unsigned int parent_index);
3720 
3727  void
3728  set_direction_flag(const bool new_direction_flag) const;
3729 
3730  template <int, int>
3731  friend class Triangulation;
3732 
3733  friend struct ::internal::TriangulationImplementation::Implementation;
3734 };
3735 
3736 
3737 
3738 /* ----- declaration of explicit specializations and general templates ----- */
3739 
3740 
3741 template <int structdim, int dim, int spacedim>
3742 template <typename OtherAccessor>
3744  const OtherAccessor &)
3745 {
3746  Assert(false,
3747  ExcMessage("You are attempting an illegal conversion between "
3748  "iterator/accessor types. The constructor you call "
3749  "only exists to make certain template constructs "
3750  "easier to write as dimension independent code but "
3751  "the conversion is not valid in the current context."));
3752 }
3753 
3754 
3755 
3756 template <int structdim, int dim, int spacedim>
3757 template <int structdim2, int dim2, int spacedim2>
3760 {
3761  Assert(false,
3762  ExcMessage("You are attempting an illegal conversion between "
3763  "iterator/accessor types. The constructor you call "
3764  "only exists to make certain template constructs "
3765  "easier to write as dimension independent code but "
3766  "the conversion is not valid in the current context."));
3767 }
3768 
3769 
3770 
3771 template <int dim, int spacedim>
3772 template <int structdim2, int dim2, int spacedim2>
3775 {
3776  Assert(false,
3777  ExcMessage("You are attempting an illegal conversion between "
3778  "iterator/accessor types. The constructor you call "
3779  "only exists to make certain template constructs "
3780  "easier to write as dimension independent code but "
3781  "the conversion is not valid in the current context."));
3782 }
3783 
3784 
3785 
3786 template <int structdim, int dim, int spacedim>
3787 template <int structdim2, int dim2, int spacedim2>
3790 {
3791  Assert(false,
3792  ExcMessage("You are attempting an illegal conversion between "
3793  "iterator/accessor types. The constructor you call "
3794  "only exists to make certain template constructs "
3795  "easier to write as dimension independent code but "
3796  "the conversion is not valid in the current context."));
3797 }
3798 
3799 
3800 
3801 template <int dim, int spacedim>
3802 template <int structdim2, int dim2, int spacedim2>
3805 {
3806  Assert(false,
3807  ExcMessage("You are attempting an illegal conversion between "
3808  "iterator/accessor types. The constructor you call "
3809  "only exists to make certain template constructs "
3810  "easier to write as dimension independent code but "
3811  "the conversion is not valid in the current context."));
3812 }
3813 
3814 
3815 #ifndef DOXYGEN
3816 
3817 template <>
3818 bool
3820 template <>
3821 bool
3823 template <>
3824 bool
3826 template <>
3827 bool
3829 template <>
3830 bool
3832 template <>
3833 bool
3835 // -------------------------------------------------------------------
3836 
3837 template <>
3838 void
3840 
3841 #endif // DOXYGEN
3842 
3844 
3845 // include more templates in debug and optimized mode
3846 #include "tria_accessor.templates.h"
3847 
3848 #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:78
const Triangulation< dim, spacedim > * tria
static ::ExceptionBase & ExcNeighborIsCoarser()
static ::ExceptionBase & ExcDereferenceInvalidObject(AccessorType arg1)
std::vector< unsigned int > vertex_indices
Definition: tria.cc:2244
bool operator!=(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
const Triangulation< 1, spacedim > * tria
unsigned int material_id
Definition: types.h:152
void set_all_manifold_ids(const types::manifold_id) const
bool operator<(const SynchronousIterators< Iterators > &a, const SynchronousIterators< Iterators > &b)
static ::ExceptionBase & ExcFacesHaveNoLevel()
bool operator==(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
static ::ExceptionBase & ExcCellNotUsed()
static ::ExceptionBase & ExcMessage(std::string arg1)
unsigned int subdomain_id
Definition: types.h:43
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:515
const Triangulation< dim, spacedim > * tria
static ::ExceptionBase & ExcSetOnlyEvenChildren(int arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1411
boost::integer_range< IncrementableType > iota_view
Definition: iota_view.h:46
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:493
#define DeclException0(Exception0)
Definition: exceptions.h:470
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
unsigned int level
Definition: grid_out.cc:4341
static ::ExceptionBase & ExcCellNotActive()
SynchronousIterators< Iterators > operator++(SynchronousIterators< Iterators > &a)
static ::ExceptionBase & ExcNoPeriodicNeighbor()
static ::ExceptionBase & ExcCantSetChildren(int arg1)
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
Definition: utilities.cc:188
Definition: cell_id.h:69
CellAccessor(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
static ::ExceptionBase & ExcCellHasNoParent()
bool point_inside(const Point< spacedim > &p) const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:362
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:152
unsigned int boundary_id
Definition: types.h:129
typename ::internal::TriaAccessorImplementation::PresentLevelType< structdim, dim >::type present_level
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcCantCompareIterators()