Reference documentation for deal.II version 9.5.0
\(\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\}}\)
Loading...
Searching...
No Matches
tria_accessor.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 1998 - 2023 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
31
32#include <boost/container/small_vector.hpp>
33
34#include <utility>
35
36
38
39// Forward declarations
40#ifndef DOXYGEN
41template <int dim, int spacedim>
43class Triangulation;
44template <typename Accessor>
45class TriaRawIterator;
46template <typename Accessor>
47class TriaIterator;
48template <typename Accessor>
50
51namespace parallel
52{
53 template <int dim, int spacedim>
55 class TriangulationBase;
56}
57
58template <int dim, int spacedim>
60class DoFHandler;
61template <int dim, int spacedim, bool lda>
62class DoFCellAccessor;
63
64
65template <int dim, int spacedim>
66class Manifold;
67
68template <int dim, int spacedim>
69class Mapping;
70#endif
71
72namespace internal
73{
74 namespace TriangulationImplementation
75 {
76 class TriaObjects;
77 struct Implementation;
78 struct ImplementationMixedMesh;
79 } // namespace TriangulationImplementation
80
81 namespace TriaAccessorImplementation
82 {
83 struct Implementation;
84
90 template <int structdim, int dim>
92 {
93 struct type
94 {
98 type() = default;
99
103 type(const int level)
104 {
106 (void)level; // removes -Wunused-parameter warning in optimized mode
107 }
108
112 operator int() const
113 {
114 return 0;
115 }
116
117 void
119 {
120 Assert(false, ExcInternalError());
121 }
122
123 void
125 {
126 Assert(false, ExcInternalError());
127 }
128 };
129 };
130
131
137 template <int dim>
138 struct PresentLevelType<dim, dim>
139 {
140 using type = int;
141 };
142 } // namespace TriaAccessorImplementation
143} // namespace internal
144template <int structdim, int dim, int spacedim>
145class TriaAccessor;
146template <int dim, int spacedim>
147class TriaAccessor<0, dim, spacedim>;
148template <int spacedim>
149class TriaAccessor<0, 1, spacedim>;
150
151// note: the file tria_accessor.templates.h is included at the end of
152// this file. this includes a lot of templates. originally, this was
153// only done in debug mode, but led to cyclic reduction problems and
154// so is now on by default.
155
156
161{
166 "The operation you are attempting can only be performed for "
167 "(cell, face, or edge) iterators that point to valid "
168 "objects. These objects need not necessarily be active, "
169 "i.e., have no children, but they need to be part of a "
170 "triangulation. (The objects pointed to by an iterator "
171 "may -- after coarsening -- also be objects that used "
172 "to be part of a triangulation, but are now no longer "
173 "used. Their memory location may have been retained "
174 "for re-use upon the next mesh refinement, but is "
175 "currently unused.)");
186 "The operation you are attempting can only be performed for "
187 "(cell, face, or edge) iterators that point to 'active' "
188 "objects. 'Active' objects are those that do not have "
189 "children (in the case of cells), or that are part of "
190 "an active cell (in the case of faces or edges). However, "
191 "the object on which you are trying the current "
192 "operation is not 'active' in this sense.");
199 "The operation you are attempting can only be performed for "
200 "(cell, face, or edge) iterators that have children, "
201 "but the object on which you are trying the current "
202 "operation does not have any.");
210 "The operation you are attempting can only be performed for "
211 "(cell, face, or edge) iterators that have a parent object, "
212 "but the object on which you are trying the current "
213 "operation does not have one -- i.e., it is on the "
214 "coarsest level of the triangulation.");
219 int,
220 << "You can only set the child index if the cell does not "
221 << "currently have children registered; or you can clear it. "
222 << "The given index was " << arg1
223 << " (-1 means: clear children).");
227 template <typename AccessorType>
229 AccessorType,
230 << "You tried to dereference an iterator for which this "
231 << "is not possible. More information on this iterator: "
232 << "index=" << arg1.index() << ", state="
233 << (arg1.state() == IteratorState::valid ?
234 "valid" :
235 (arg1.state() == IteratorState::past_the_end ?
236 "past_the_end" :
237 "invalid")));
242 "Iterators can only be compared if they point to the same "
243 "triangulation, or if neither of them are associated "
244 "with a triangulation.");
245 // TODO: Write documentation!
250 // TODO: Write documentation!
270 // TODO: Write documentation!
276 int,
277 << "You can only set the child index of an even numbered child."
278 << "The number of the child given was " << arg1 << '.');
279} // namespace TriaAccessorExceptions
280
281
307template <int structdim, int dim, int spacedim = dim>
309{
310public:
317 static constexpr unsigned int space_dimension = spacedim;
318
324 static constexpr unsigned int dimension = dim;
325
331 static const unsigned int structure_dimension = structdim;
332
342 void
343 operator=(const TriaAccessorBase *) = delete;
344
345protected:
351 using AccessorData = void;
352
357 const int level = -1,
358 const int index = -1,
359 const AccessorData * = nullptr);
360
365
373 void
375
381
392 bool
393 operator<(const TriaAccessorBase &other) const;
394
395protected:
399 bool
401
405 bool
407
421 void
423
431 void
441 objects() const;
442
443public:
449 using LocalData = void *;
450
474 int
475 level() const;
476
503 int
504 index() const;
505
511 state() const;
512
519
523protected:
528 typename ::internal::TriaAccessorImplementation::
529 PresentLevelType<structdim, dim>::type present_level;
530
536
541
542private:
543 template <typename Accessor>
544 friend class TriaRawIterator;
545 template <typename Accessor>
546 friend class TriaIterator;
547 template <typename Accessor>
548 friend class TriaActiveIterator;
549};
550
551
552
573template <int structdim, int dim, int spacedim = dim>
575{
576public:
583 static constexpr unsigned int space_dimension = spacedim;
584
590 static constexpr unsigned int dimension = dim;
591
597 static const unsigned int structure_dimension = structdim;
598
604 using AccessorData = void;
605
613 InvalidAccessor(const void * parent = nullptr,
614 const int level = -1,
615 const int index = -1,
616 const AccessorData *local_data = nullptr);
617
626
631 template <typename OtherAccessor>
632 InvalidAccessor(const OtherAccessor &);
633
637 void
639
643 bool
645 bool
647
651 void
652 operator++() const;
653 void
654 operator--() const;
655
661 state();
662
663
668 static int
669 level();
670
675 static int
676 index();
677
682 bool
683 used() const;
684
689 bool
691
696 manifold_id() const;
697
701 unsigned int
702 user_index() const;
703
707 void
708 set_user_index(const unsigned int p) const;
709
713 void
715
720 vertex(const unsigned int i) const;
721
726 void *
727 line(const unsigned int i) const;
728
733 void *
734 quad(const unsigned int i) const;
735};
736
737
738
756template <int structdim, int dim, int spacedim>
757class TriaAccessor : public TriaAccessorBase<structdim, dim, spacedim>
758{
759public:
765
770 const int level = -1,
771 const int index = -1,
772 const AccessorData * local_data = nullptr);
773
778 TriaAccessor(const TriaAccessor &) = default;
779
783 TriaAccessor(TriaAccessor &&) = default; // NOLINT
784
797 template <int structdim2, int dim2, int spacedim2>
799
804 template <int structdim2, int dim2, int spacedim2>
806
817 operator=(const TriaAccessor &) = delete;
818
823 operator=(TriaAccessor &&) = default; // NOLINT
824
828 ~TriaAccessor() = default;
829
836 bool
837 used() const;
838
851 vertex_iterator(const unsigned int i) const;
852
868 unsigned int
869 vertex_index(const unsigned int i) const;
870
909 vertex(const unsigned int i) const;
910
914 typename ::internal::TriangulationImplementation::
915 Iterators<dim, spacedim>::line_iterator
916 line(const unsigned int i) const;
917
924 unsigned int
925 line_index(const unsigned int i) const;
926
930 typename ::internal::TriangulationImplementation::
931 Iterators<dim, spacedim>::quad_iterator
932 quad(const unsigned int i) const;
933
940 unsigned int
941 quad_index(const unsigned int i) const;
957 unsigned char
958 combined_face_orientation(const unsigned int face) const;
959
971 bool
972 face_orientation(const unsigned int face) const;
973
983 bool
984 face_flip(const unsigned int face) const;
985
995 bool
996 face_rotation(const unsigned int face) const;
997
1010 bool
1011 line_orientation(const unsigned int line) const;
1026 bool
1028
1033 unsigned int
1034 n_children() const;
1035
1040 unsigned int
1042
1056 unsigned int
1058
1072 unsigned int
1074
1079 child(const unsigned int i) const;
1080
1085 unsigned int
1088
1098 isotropic_child(const unsigned int i) const;
1099
1105
1111 int
1112 child_index(const unsigned int i) const;
1113
1119 int
1120 isotropic_child_index(const unsigned int i) const;
1144
1174 void
1176
1207 void
1209
1217 bool
1219
1231
1254
1272 void
1274
1288 void
1290
1307 bool
1309
1315 void
1317
1323 void
1325
1331 void
1333
1339 void
1341
1347 void
1349
1361 void
1362 set_user_pointer(void *p) const;
1363
1369 void
1371
1387 void *
1389
1411 void
1413
1420 void
1422
1432 void
1433 set_user_index(const unsigned int p) const;
1434
1440 void
1442
1454 unsigned int
1455 user_index() const;
1456
1474 void
1475 recursively_set_user_index(const unsigned int p) const;
1476
1485 void
1521 double
1522 diameter() const;
1523
1550 std::pair<Point<spacedim>, double>
1552
1562 bounding_box() const;
1563
1573 double
1574 extent_in_direction(const unsigned int axis) const;
1575
1579 double
1581
1596 intermediate_point(const Point<structdim> &coordinates) const;
1597
1622
1658 center(const bool respect_manifold = false,
1659 const bool interpolate_from_surrounding = false) const;
1660
1679 barycenter() const;
1680
1706 double
1707 measure() const;
1708
1723 bool
1726
1732
1736 unsigned int
1737 n_vertices() const;
1738
1742 unsigned int
1743 n_lines() const;
1744
1754 unsigned int
1755 n_faces() const;
1756
1763
1770
1779
1784private:
1789 void
1791
1799 void
1801 const std::initializer_list<int> &new_indices) const;
1802
1806 void
1808 const std::initializer_list<unsigned int> &new_indices) const;
1809
1817 void
1818 set_line_orientation(const unsigned int line, const bool orientation) const;
1819
1827 void
1828 set_combined_face_orientation(const unsigned int face,
1829 const unsigned char combined_orientation) const;
1830
1834 void
1836
1840 void
1842
1851 void
1853
1861 void
1863
1870 void
1871 set_children(const unsigned int i, const int index) const;
1872
1877 void
1879
1880private:
1881 friend class Triangulation<dim, spacedim>;
1882
1883 friend struct ::internal::TriangulationImplementation::Implementation;
1884 friend struct ::internal::TriangulationImplementation::
1885 ImplementationMixedMesh;
1886 friend struct ::internal::TriaAccessorImplementation::Implementation;
1887};
1888
1889
1890
1909template <int dim, int spacedim>
1910class TriaAccessor<0, dim, spacedim>
1911{
1912public:
1918 static constexpr unsigned int space_dimension = spacedim;
1919
1925 static constexpr unsigned int dimension = dim;
1926
1932 static const unsigned int structure_dimension = 0;
1933
1937 using AccessorData = void;
1938
1944 const unsigned int vertex_index);
1945
1952 const int level = 0,
1953 const int index = 0,
1954 const AccessorData * = nullptr);
1955
1959 template <int structdim2, int dim2, int spacedim2>
1961
1965 template <int structdim2, int dim2, int spacedim2>
1967
1972 state() const;
1973
1978 static int
1980
1985 int
1986 index() const;
1987
1994
2004 void
2006
2010 void
2015 bool
2016 operator==(const TriaAccessor &) const;
2017
2021 bool
2022 operator!=(const TriaAccessor &) const;
2023
2051 unsigned int
2052 vertex_index(const unsigned int i = 0) const;
2053
2060 vertex(const unsigned int i = 0) const;
2061
2066 typename ::internal::TriangulationImplementation::
2067 Iterators<dim, spacedim>::line_iterator static line(const unsigned int);
2068
2072 static unsigned int
2073 line_index(const unsigned int i);
2074
2078 static typename ::internal::TriangulationImplementation::
2079 Iterators<dim, spacedim>::quad_iterator
2080 quad(const unsigned int i);
2081
2085 static unsigned int
2086 quad_index(const unsigned int i);
2087
2103 double
2104 diameter() const;
2105
2113 double
2114 extent_in_direction(const unsigned int axis) const;
2115
2124 center(const bool respect_manifold = false,
2125 const bool interpolate_from_surrounding = false) const;
2126
2135 double
2136 measure() const;
2151 static unsigned char
2152 combined_face_orientation(const unsigned int face);
2153
2157 static bool
2158 face_orientation(const unsigned int face);
2159
2163 static bool
2164 face_flip(const unsigned int face);
2165
2169 static bool
2170 face_rotation(const unsigned int face);
2171
2175 static bool
2176 line_orientation(const unsigned int line);
2177
2192 static bool
2194
2199 static unsigned int
2201
2206 static unsigned int
2208
2213 static unsigned int
2215
2216
2220 static unsigned int
2222
2226 static unsigned int
2228
2233 child(const unsigned int);
2234
2239 isotropic_child(const unsigned int);
2240
2244 static RefinementCase<0>
2246
2250 static int
2251 child_index(const unsigned int i);
2252
2256 static int
2257 isotropic_child_index(const unsigned int i);
2265 bool
2266 used() const;
2267
2268protected:
2276 void
2278
2287 bool
2288 operator<(const TriaAccessor &other) const;
2289
2294
2299
2300private:
2301 template <typename Accessor>
2302 friend class TriaRawIterator;
2303 template <typename Accessor>
2304 friend class TriaIterator;
2305 template <typename Accessor>
2307};
2308
2309
2310
2327template <int spacedim>
2328class TriaAccessor<0, 1, spacedim>
2329{
2330public:
2336 static constexpr unsigned int space_dimension = spacedim;
2337
2343 static constexpr unsigned int dimension = 1;
2344
2350 static const unsigned int structure_dimension = 0;
2351
2355 using AccessorData = void;
2356
2362 {
2374 right_vertex
2376
2389 const VertexKind vertex_kind,
2390 const unsigned int vertex_index);
2391
2398 const int = 0,
2399 const int = 0,
2400 const AccessorData * = nullptr);
2401
2405 template <int structdim2, int dim2, int spacedim2>
2407
2411 template <int structdim2, int dim2, int spacedim2>
2413
2418 void
2420
2426 void
2428
2436
2441 static int
2443
2448 int
2449 index() const;
2450
2457
2468 void
2469 operator++() const;
2470
2475 void
2476 operator--() const;
2480 bool
2481 operator==(const TriaAccessor &) const;
2482
2486 bool
2487 operator!=(const TriaAccessor &) const;
2488
2497 bool
2498 operator<(const TriaAccessor &other) const;
2499
2526 unsigned int
2527 vertex_index(const unsigned int i = 0) const;
2528
2535 vertex(const unsigned int i = 0) const;
2536
2542 center() const;
2543
2548 typename ::internal::TriangulationImplementation::
2549 Iterators<1, spacedim>::line_iterator static line(const unsigned int);
2550
2557 static unsigned int
2558 line_index(const unsigned int i);
2559
2563 static typename ::internal::TriangulationImplementation::
2564 Iterators<1, spacedim>::quad_iterator
2565 quad(const unsigned int i);
2566
2573 static unsigned int
2574 quad_index(const unsigned int i);
2575
2585 bool
2587
2604
2608 const Manifold<1, spacedim> &
2610
2619
2620
2632 bool
2633 user_flag_set() const;
2634
2640 void
2641 set_user_flag() const;
2642
2648 void
2649 clear_user_flag() const;
2650
2656 void
2658
2664 void
2666
2672 void
2673 clear_user_data() const;
2674
2686 void
2687 set_user_pointer(void *p) const;
2688
2694 void
2695 clear_user_pointer() const;
2696
2712 void *
2713 user_pointer() const;
2714
2736 void
2737 recursively_set_user_pointer(void *p) const;
2738
2745 void
2747
2757 void
2758 set_user_index(const unsigned int p) const;
2759
2765 void
2766 clear_user_index() const;
2767
2779 unsigned int
2780 user_index() const;
2781
2799 void
2800 recursively_set_user_index(const unsigned int p) const;
2801
2810 void
2826 static unsigned char
2827 combined_face_orientation(const unsigned int face);
2828
2832 static bool
2833 face_orientation(const unsigned int face);
2834
2838 static bool
2839 face_flip(const unsigned int face);
2840
2844 static bool
2845 face_rotation(const unsigned int face);
2846
2850 static bool
2851 line_orientation(const unsigned int line);
2852
2867 static bool
2869
2874 static unsigned int
2876
2881 static unsigned int
2883
2888 static unsigned int
2890
2891
2895 static unsigned int
2897
2901 static unsigned int
2903
2908 child(const unsigned int);
2909
2914 isotropic_child(const unsigned int);
2915
2919 static RefinementCase<0>
2921
2925 static int
2926 child_index(const unsigned int i);
2927
2931 static int
2932 isotropic_child_index(const unsigned int i);
2965 void
2967
2974 void
2976
2988 void
2990
3002 void
3011 bool
3012 used() const;
3013
3019
3023 unsigned int
3024 n_vertices() const;
3025
3029 unsigned int
3030 n_lines() const;
3031
3038
3045
3046protected:
3051
3057
3062};
3063
3064
3065
3081template <int dim, int spacedim = dim>
3082class CellAccessor : public TriaAccessor<dim, dim, spacedim>
3083{
3084public:
3089
3094
3106 const int level = -1,
3107 const int index = -1,
3108 const AccessorData * local_data = nullptr);
3109
3114
3127 template <int structdim2, int dim2, int spacedim2>
3129
3134 template <int structdim2, int dim2, int spacedim2>
3136
3141
3145 // NOLINTNEXTLINE OSX does not compile with noexcept
3147
3151 ~CellAccessor() = default;
3152
3164
3168 // NOLINTNEXTLINE OSX does not compile with noexcept
3171
3195 as_dof_handler_iterator(const DoFHandler<dim, spacedim> &dof_handler) const;
3196
3213 child(const unsigned int i) const;
3214
3218 boost::container::small_vector<TriaIterator<CellAccessor<dim, spacedim>>,
3221
3225 TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>
3226 face(const unsigned int i) const;
3227
3232 unsigned int
3235
3239 boost::container::small_vector<
3240 TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>,
3243
3253 unsigned int
3254 face_index(const unsigned int i) const;
3255
3304 neighbor_child_on_subface(const unsigned int face_no,
3305 const unsigned int subface_no) const;
3306
3332 neighbor(const unsigned int face_no) const;
3333
3341 int
3342 neighbor_index(const unsigned int face_no) const;
3343
3351 int
3352 neighbor_level(const unsigned int face_no) const;
3353
3365 unsigned int
3366 neighbor_of_neighbor(const unsigned int face_no) const;
3367
3378 bool
3379 neighbor_is_coarser(const unsigned int face_no) const;
3380
3395 std::pair<unsigned int, unsigned int>
3396 neighbor_of_coarser_neighbor(const unsigned int neighbor) const;
3397
3404 unsigned int
3405 neighbor_face_no(const unsigned int neighbor) const;
3406
3410 static bool
3412
3426 bool
3427 has_periodic_neighbor(const unsigned int i) const;
3428
3446 periodic_neighbor(const unsigned int i) const;
3447
3456 neighbor_or_periodic_neighbor(const unsigned int i) const;
3457
3473 periodic_neighbor_child_on_subface(const unsigned int face_no,
3474 const unsigned int subface_no) const;
3475
3486 std::pair<unsigned int, unsigned int>
3487 periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const;
3488
3494 int
3495 periodic_neighbor_index(const unsigned int i) const;
3496
3502 int
3503 periodic_neighbor_level(const unsigned int i) const;
3504
3519 unsigned int
3520 periodic_neighbor_of_periodic_neighbor(const unsigned int i) const;
3521
3527 unsigned int
3528 periodic_neighbor_face_no(const unsigned int i) const;
3529
3536 bool
3537 periodic_neighbor_is_coarser(const unsigned int i) const;
3538
3555 bool
3556 at_boundary(const unsigned int i) const;
3557
3566 bool
3567 at_boundary() const;
3568
3576 bool
3577 has_boundary_lines() const;
3605
3623 void
3626
3630 void
3632
3640 bool
3642 const unsigned int face_no,
3643 const RefinementCase<dim - 1> &face_refinement_case =
3645
3651 bool
3652 flag_for_line_refinement(const unsigned int line_no) const;
3653
3663 subface_case(const unsigned int face_no) const;
3664
3668 bool
3670
3675 void
3677
3681 void
3706 material_id() const;
3707
3719 void
3720 set_material_id(const types::material_id new_material_id) const;
3721
3730 void
3731 recursively_set_material_id(const types::material_id new_material_id) const;
3759
3775 void
3776 set_subdomain_id(const types::subdomain_id new_subdomain_id) const;
3777
3788
3793 void
3795 const types::subdomain_id new_level_subdomain_id) const;
3796
3797
3813 void
3815 const types::subdomain_id new_subdomain_id) const;
3839
3849
3863 bool
3864 direction_flag() const;
3865
3891 unsigned int
3893
3901 int
3902 parent_index() const;
3903
3910 parent() const;
3911
3931 bool
3932 is_active() const;
3933
3953 bool
3955
3960 bool
3962
3996 bool
3997 is_ghost() const;
3998
4004 bool
4006
4033 bool
4035
4042 bool
4044
4058 bool
4060
4069 void
4070 set_neighbor(const unsigned int i,
4071 const TriaIterator<CellAccessor<dim, spacedim>> &pointer) const;
4072
4086 CellId
4087 id() const;
4088
4089 using TriaAccessor<dim, dim, spacedim>::diameter;
4090
4094 double
4095 diameter(const Mapping<dim, spacedim> &mapping) const;
4096
4114
4115protected:
4131 unsigned int
4132 neighbor_of_neighbor_internal(const unsigned int neighbor) const;
4133
4139 template <int dim_, int spacedim_>
4140 bool
4141 point_inside_codim(const Point<spacedim_> &p) const;
4142
4143
4144
4145private:
4150 void
4151 set_active_cell_index(const unsigned int active_cell_index) const;
4152
4156 void
4158
4162 void
4164
4168 void
4169 set_parent(const unsigned int parent_index);
4170
4177 void
4178 set_direction_flag(const bool new_direction_flag) const;
4179
4180 friend class Triangulation<dim, spacedim>;
4181
4182 friend class parallel::TriangulationBase<dim, spacedim>;
4183
4184 friend struct ::internal::TriangulationImplementation::Implementation;
4185 friend struct ::internal::TriangulationImplementation::
4186 ImplementationMixedMesh;
4187};
4188
4189
4190
4191/* ----- declaration of explicit specializations and general templates ----- */
4192
4193
4194template <int structdim, int dim, int spacedim>
4195template <typename OtherAccessor>
4197 const OtherAccessor &)
4198{
4199 Assert(false,
4200 ExcMessage("You are attempting an illegal conversion between "
4201 "iterator/accessor types. The constructor you call "
4202 "only exists to make certain template constructs "
4203 "easier to write as dimension independent code but "
4204 "the conversion is not valid in the current context."));
4205}
4206
4207
4208
4209template <int structdim, int dim, int spacedim>
4210template <int structdim2, int dim2, int spacedim2>
4213{
4214 Assert(false,
4215 ExcMessage("You are attempting an illegal conversion between "
4216 "iterator/accessor types. The constructor you call "
4217 "only exists to make certain template constructs "
4218 "easier to write as dimension independent code but "
4219 "the conversion is not valid in the current context."));
4220}
4221
4222
4223
4224template <int dim, int spacedim>
4225template <int structdim2, int dim2, int spacedim2>
4228{
4229 Assert(false,
4230 ExcMessage("You are attempting an illegal conversion between "
4231 "iterator/accessor types. The constructor you call "
4232 "only exists to make certain template constructs "
4233 "easier to write as dimension independent code but "
4234 "the conversion is not valid in the current context."));
4235}
4236
4237
4238
4239template <int structdim, int dim, int spacedim>
4240template <int structdim2, int dim2, int spacedim2>
4243{
4244 Assert(false,
4245 ExcMessage("You are attempting an illegal conversion between "
4246 "iterator/accessor types. The constructor you call "
4247 "only exists to make certain template constructs "
4248 "easier to write as dimension independent code but "
4249 "the conversion is not valid in the current context."));
4250}
4251
4252
4253
4254template <int dim, int spacedim>
4255template <int structdim2, int dim2, int spacedim2>
4258{
4259 Assert(false,
4260 ExcMessage("You are attempting an illegal conversion between "
4261 "iterator/accessor types. The constructor you call "
4262 "only exists to make certain template constructs "
4263 "easier to write as dimension independent code but "
4264 "the conversion is not valid in the current context."));
4265}
4266
4267
4268#ifndef DOXYGEN
4269
4270template <>
4271bool
4273template <>
4274bool
4276template <>
4277bool
4279template <>
4280bool
4282template <>
4283bool
4285template <>
4286bool
4288// -------------------------------------------------------------------
4289
4290template <>
4291void
4293
4294#endif // DOXYGEN
4295
4297
4298// include more templates in debug and optimized mode
4299#include <deal.II/grid/tria_accessor.templates.h>
4300
4301#endif
void recursively_set_subdomain_id(const types::subdomain_id new_subdomain_id) const
CellAccessor(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
TriaIterator< CellAccessor< dim, spacedim > > parent() const
unsigned int neighbor_face_no(const unsigned int neighbor) const
void set_active_cell_index(const unsigned int active_cell_index) const
unsigned int neighbor_of_neighbor_internal(const unsigned int neighbor) const
TriaIterator< CellAccessor< dim, spacedim > > periodic_neighbor(const unsigned int i) const
types::global_cell_index global_active_cell_index() const
void set_neighbor(const unsigned int i, const TriaIterator< CellAccessor< dim, spacedim > > &pointer) const
TriaIterator< CellAccessor< dim, spacedim > > neighbor(const unsigned int face_no) const
bool is_artificial_on_level() const
void set_direction_flag(const bool new_direction_flag) const
void recursively_set_material_id(const types::material_id new_material_id) const
void set_level_subdomain_id(const types::subdomain_id new_level_subdomain_id) const
types::subdomain_id level_subdomain_id() const
RefinementCase< dim > refine_flag_set() const
bool flag_for_face_refinement(const unsigned int face_no, const RefinementCase< dim - 1 > &face_refinement_case=RefinementCase< dim - 1 >::isotropic_refinement) const
unsigned int face_index(const unsigned int i) const
double diameter(const Mapping< dim, spacedim > &mapping) const
TriaActiveIterator< DoFCellAccessor< dim, spacedim, false > > as_dof_handler_iterator(const DoFHandler< dim, spacedim > &dof_handler) const
CellAccessor< dim, spacedim > & operator=(CellAccessor< dim, spacedim > &&)=default
bool is_active() const
bool is_ghost() const
TriaIterator< CellAccessor< dim, spacedim > > neighbor_child_on_subface(const unsigned int face_no, const unsigned int subface_no) const
void set_subdomain_id(const types::subdomain_id new_subdomain_id) const
bool neighbor_is_coarser(const unsigned int face_no) const
TriaIterator< CellAccessor< dim, spacedim > > child(const unsigned int i) const
void set_global_level_cell_index(const types::global_cell_index index) const
bool has_periodic_neighbor(const unsigned int i) const
int periodic_neighbor_level(const unsigned int i) const
std::pair< unsigned int, unsigned int > neighbor_of_coarser_neighbor(const unsigned int neighbor) const
~CellAccessor()=default
CellAccessor(const CellAccessor< dim, spacedim > &)=default
void set_coarsen_flag() const
TriaIterator< TriaAccessor< dim - 1, dim, spacedim > > face(const unsigned int i) const
CellAccessor< dim, spacedim > & operator=(const CellAccessor< dim, spacedim > &)=delete
boost::container::small_vector< TriaIterator< TriaAccessor< dim - 1, dim, spacedim > >, GeometryInfo< dim >::faces_per_cell > face_iterators() const
unsigned int neighbor_of_neighbor(const unsigned int face_no) const
unsigned int face_iterator_to_index(const TriaIterator< TriaAccessor< dim - 1, dim, spacedim > > &face) const
void set_material_id(const types::material_id new_material_id) const
bool is_locally_owned() const
void set_refine_flag(const RefinementCase< dim > ref_case=RefinementCase< dim >::isotropic_refinement) const
CellAccessor(CellAccessor< dim, spacedim > &&)=default
bool point_inside_codim(const Point< spacedim_ > &p) const
typename TriaAccessor< dim, dim, spacedim >::AccessorData AccessorData
bool is_locally_owned_on_level() const
bool has_boundary_lines() const
TriaIterator< CellAccessor< dim, spacedim > > periodic_neighbor_child_on_subface(const unsigned int face_no, const unsigned int subface_no) const
int neighbor_level(const unsigned int face_no) const
int periodic_neighbor_index(const unsigned int i) const
bool periodic_neighbor_is_coarser(const unsigned int i) const
void set_global_active_cell_index(const types::global_cell_index index) const
void clear_coarsen_flag() const
boost::container::small_vector< TriaIterator< CellAccessor< dim, spacedim > >, GeometryInfo< dim >::max_children_per_cell > child_iterators() const
void set_parent(const unsigned int parent_index)
std::pair< unsigned int, unsigned int > periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const
CellAccessor(const TriaAccessor< dim, dim, spacedim > &cell_accessor)
bool at_boundary() const
unsigned int active_cell_index() const
static bool is_level_cell()
void clear_refine_flag() const
int neighbor_index(const unsigned int face_no) const
bool point_inside(const Point< spacedim > &p) const
types::subdomain_id subdomain_id() const
bool direction_flag() const
types::material_id material_id() const
bool coarsen_flag_set() const
types::global_cell_index global_level_cell_index() const
bool flag_for_line_refinement(const unsigned int line_no) const
CellId id() const
::internal::SubfaceCase< dim > subface_case(const unsigned int face_no) const
bool is_artificial() const
bool is_ghost_on_level() const
TriaIterator< CellAccessor< dim, spacedim > > neighbor_or_periodic_neighbor(const unsigned int i) const
int parent_index() const
unsigned int periodic_neighbor_of_periodic_neighbor(const unsigned int i) const
unsigned int periodic_neighbor_face_no(const unsigned int i) const
void set_manifold_id(const types::manifold_id) const
InvalidAccessor(const InvalidAccessor &)
static constexpr unsigned int space_dimension
static int level()
bool operator==(const InvalidAccessor &) const
Point< spacedim > & vertex(const unsigned int i) const
void * quad(const unsigned int i) const
void operator++() const
void * line(const unsigned int i) const
static IteratorState::IteratorStates state()
bool operator!=(const InvalidAccessor &) const
bool has_children() const
unsigned int user_index() const
static const unsigned int structure_dimension
void copy_from(const InvalidAccessor &)
void set_user_index(const unsigned int p) const
static int index()
types::manifold_id manifold_id() const
InvalidAccessor(const void *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
static constexpr unsigned int dimension
void operator--() const
bool used() const
InvalidAccessor(const OtherAccessor &)
Abstract base class for mapping classes.
Definition mapping.h:317
Definition point.h:112
bool operator!=(const TriaAccessorBase &) const
static constexpr unsigned int space_dimension
static constexpr unsigned int dimension
TriaAccessorBase(const TriaAccessorBase &)
void operator=(const TriaAccessorBase *)=delete
static const unsigned int structure_dimension
TriaAccessorBase(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *=nullptr)
void copy_from(const TriaAccessorBase &)
const Triangulation< dim, spacedim > & get_triangulation() const
IteratorState::IteratorStates state() const
int index() const
bool operator<(const TriaAccessorBase &other) const
typename::internal::TriaAccessorImplementation::PresentLevelType< structdim, dim >::type present_level
::internal::TriangulationImplementation::TriaObjects & objects() const
TriaAccessorBase & operator=(const TriaAccessorBase &)
int level() const
const Triangulation< dim, spacedim > * tria
bool operator==(const TriaAccessorBase &) const
static TriaIterator< TriaAccessor< 0, 1, spacedim > > isotropic_child(const unsigned int)
Return an invalid object.
static unsigned int line_index(const unsigned int i)
std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices() const
static bool face_orientation(const unsigned int face)
Always return false.
static unsigned int number_of_children()
static unsigned int child_iterator_to_index(const TriaIterator< TriaAccessor< 0, 1, spacedim > > &)
Return an invalid unsigned integer.
static unsigned int max_refinement_depth()
static unsigned int quad_index(const unsigned int i)
static typename::internal::TriangulationImplementation::Iterators< 1, spacedim >::quad_iterator quad(const unsigned int i)
const Manifold< 1, spacedim > & get_manifold() const
TriaAccessor(const InvalidAccessor< structdim2, dim2, spacedim2 > &)
static bool face_flip(const unsigned int face)
Always return false.
TriaAccessor(const Triangulation< 1, spacedim > *tria, const VertexKind vertex_kind, const unsigned int vertex_index)
static TriaIterator< TriaAccessor< 0, 1, spacedim > > child(const unsigned int)
Return an invalid object.
unsigned int n_lines() const
const Triangulation< 1, spacedim > & get_triangulation() const
bool operator!=(const TriaAccessor &) const
TriaAccessor(const Triangulation< 1, spacedim > *tria=nullptr, const int=0, const int=0, const AccessorData *=nullptr)
static unsigned int n_active_descendants()
static int isotropic_child_index(const unsigned int i)
Returns -1.
types::boundary_id boundary_id() const
static bool line_orientation(const unsigned int line)
Always return false.
static typename::internal::TriangulationImplementation::Iterators< 1, spacedim >::line_iterator line(const unsigned int)
static RefinementCase< 0 > refinement_case()
ReferenceCell reference_cell() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > line_indices() const
unsigned int vertex_index(const unsigned int i=0) const
void copy_from(const TriaAccessor &)
void set_all_boundary_ids(const types::boundary_id) const
static IteratorState::IteratorStates state()
static bool face_rotation(const unsigned int face)
Always return false.
Point< spacedim > center() const
static int child_index(const unsigned int i)
Returns -1.
types::manifold_id manifold_id() const
static unsigned int n_children()
TriaAccessor(const TriaAccessor< structdim2, dim2, spacedim2 > &)
static unsigned char combined_face_orientation(const unsigned int face)
Always return 0.
void set_boundary_id(const types::boundary_id) const
Point< spacedim > & vertex(const unsigned int i=0) const
unsigned int n_vertices() const
void set_manifold_id(const types::manifold_id)
void copy_from(const TriaAccessorBase< 0, 1, spacedim > &)
bool operator==(const TriaAccessor &) const
const Triangulation< 1, spacedim > * tria
double extent_in_direction(const unsigned int axis) const
TriaAccessor(const Triangulation< dim, spacedim > *tria, const unsigned int vertex_index)
Point< spacedim > & vertex(const unsigned int i=0) const
static typename::internal::TriangulationImplementation::Iterators< dim, spacedim >::line_iterator line(const unsigned int)
bool operator!=(const TriaAccessor &) const
static RefinementCase< 0 > refinement_case()
void copy_from(const TriaAccessor &)
static TriaIterator< TriaAccessor< 0, dim, spacedim > > child(const unsigned int)
Return an invalid object.
TriaAccessor(const InvalidAccessor< structdim2, dim2, spacedim2 > &)
static int child_index(const unsigned int i)
Returns -1.
bool operator==(const TriaAccessor &) const
static unsigned int line_index(const unsigned int i)
static TriaIterator< TriaAccessor< 0, dim, spacedim > > isotropic_child(const unsigned int)
Return an invalid object.
static typename::internal::TriangulationImplementation::Iterators< dim, spacedim >::quad_iterator quad(const unsigned int i)
unsigned int vertex_index(const unsigned int i=0) const
const Triangulation< dim, spacedim > * tria
static int isotropic_child_index(const unsigned int i)
Returns -1.
TriaAccessor(const Triangulation< dim, spacedim > *tria=nullptr, const int level=0, const int index=0, const AccessorData *=nullptr)
static unsigned char combined_face_orientation(const unsigned int face)
Always return 0.
static unsigned int n_children()
const Triangulation< dim, spacedim > & get_triangulation() const
static unsigned int child_iterator_to_index(const TriaIterator< TriaAccessor< 0, dim, spacedim > > &)
Return an invalid unsigned integer.
static bool face_flip(const unsigned int face)
Always return false.
static bool line_orientation(const unsigned int line)
Always return false.
static unsigned int max_refinement_depth()
static unsigned int quad_index(const unsigned int i)
static bool face_rotation(const unsigned int face)
Always return false.
static bool face_orientation(const unsigned int face)
Always return false.
static unsigned int n_active_descendants()
TriaAccessor(const TriaAccessor< structdim2, dim2, spacedim2 > &)
Point< spacedim > center(const bool respect_manifold=false, const bool interpolate_from_surrounding=false) const
IteratorState::IteratorStates state() const
static unsigned int number_of_children()
std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices() const
void clear_children() const
bool line_orientation(const unsigned int line) const
void set_boundary_id_internal(const types::boundary_id id) const
void set_user_index(const unsigned int p) const
TriaAccessor(TriaAccessor &&)=default
void clear_user_pointer() const
unsigned int n_active_descendants() const
unsigned char combined_face_orientation(const unsigned int face) const
bool face_rotation(const unsigned int face) const
void recursively_set_user_index(const unsigned int p) const
TriaIterator< TriaAccessor< 0, dim, spacedim > > vertex_iterator(const unsigned int i) const
void clear_user_data() const
TriaAccessor(const TriaAccessor< structdim2, dim2, spacedim2 > &)
Point< structdim > real_to_unit_cell_affine_approximation(const Point< spacedim > &point) const
void recursively_clear_user_index() const
unsigned int line_index(const unsigned int i) const
bool face_orientation(const unsigned int face) const
const Manifold< dim, spacedim > & get_manifold() const
void recursively_set_user_pointer(void *p) const
unsigned int n_lines() const
double extent_in_direction(const unsigned int axis) const
Point< spacedim > intermediate_point(const Point< structdim > &coordinates) const
unsigned int n_vertices() const
bool has_children() const
void recursively_clear_user_flag() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > line_indices() const
typename TriaAccessorBase< structdim, dim, spacedim >::AccessorData AccessorData
Point< spacedim > barycenter() const
BoundingBox< spacedim > bounding_box() const
void clear_user_flag() const
TriaIterator< TriaAccessor< structdim, dim, spacedim > > child(const unsigned int i) const
unsigned int n_children() const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > face_indices() const
void recursively_set_user_flag() const
void set_boundary_id(const types::boundary_id) const
bool is_translation_of(const TriaIterator< TriaAccessor< structdim, dim, spacedim > > &o) const
bool user_flag_set() const
void set_used_flag() const
types::manifold_id manifold_id() const
void set_user_flag() const
TriaAccessor(const TriaAccessor &)=default
void clear_used_flag() const
std::pair< Point< spacedim >, double > enclosing_ball() const
unsigned int vertex_index(const unsigned int i) const
void * user_pointer() const
int isotropic_child_index(const unsigned int i) const
~TriaAccessor()=default
void set_refinement_case(const RefinementCase< structdim > &ref_case) const
void clear_user_index() const
double minimum_vertex_distance() const
double measure() const
void set_all_boundary_ids(const types::boundary_id) const
void set_bounding_object_indices(const std::initializer_list< int > &new_indices) const
unsigned int number_of_children() const
unsigned int max_refinement_depth() const
Point< spacedim > & vertex(const unsigned int i) const
void set_line_orientation(const unsigned int line, const bool orientation) const
unsigned int quad_index(const unsigned int i) const
TriaAccessor & operator=(const TriaAccessor &)=delete
unsigned int user_index() const
int child_index(const unsigned int i) const
void set_user_pointer(void *p) const
void recursively_clear_user_pointer() const
ReferenceCell reference_cell() const
RefinementCase< structdim > refinement_case() const
typename::internal::TriangulationImplementation::Iterators< dim, spacedim >::line_iterator line(const unsigned int i) const
bool face_flip(const unsigned int face) const
TriaAccessor(const Triangulation< dim, spacedim > *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
unsigned int child_iterator_to_index(const TriaIterator< TriaAccessor< structdim, dim, spacedim > > &child) const
TriaAccessor(const InvalidAccessor< structdim2, dim2, spacedim2 > &)
bool used() const
void set_children(const unsigned int i, const int index) const
void clear_refinement_case() const
double diameter() const
bool at_boundary() const
types::boundary_id boundary_id() const
unsigned int n_faces() const
TriaAccessor & operator=(TriaAccessor &&)=default
TriaIterator< TriaAccessor< structdim, dim, spacedim > > isotropic_child(const unsigned int i) const
typename::internal::TriangulationImplementation::Iterators< dim, spacedim >::quad_iterator quad(const unsigned int i) const
void set_combined_face_orientation(const unsigned int face, const unsigned char combined_orientation) const
#define DEAL_II_DEPRECATED
Definition config.h:172
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_CXX20_REQUIRES(condition)
Definition config.h:160
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
Point< 3 > center
unsigned int level
Definition grid_out.cc:4618
static ::ExceptionBase & ExcCantCompareIterators()
#define DeclException0(Exception0)
Definition exceptions.h:465
static ::ExceptionBase & ExcCellHasNoParent()
static ::ExceptionBase & ExcCellFlaggedForRefinement()
static ::ExceptionBase & ExcCellHasNoChildren()
static ::ExceptionBase & ExcCellNotUsed()
static ::ExceptionBase & ExcNoPeriodicNeighbor()
static ::ExceptionBase & ExcCantSetChildren(int arg1)
#define Assert(cond, exc)
static ::ExceptionBase & ExcNeighborIsNotCoarser()
static ::ExceptionBase & ExcSetOnlyEvenChildren(int arg1)
static ::ExceptionBase & ExcNeighborIsCoarser()
#define DeclExceptionMsg(Exception, defaulttext)
Definition exceptions.h:488
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcFacesHaveNoLevel()
static ::ExceptionBase & ExcDereferenceInvalidObject(AccessorType arg1)
static ::ExceptionBase & ExcCellFlaggedForCoarsening()
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:510
static ::ExceptionBase & ExcRefineCellNotActive()
static ::ExceptionBase & ExcMessage(std::string arg1)
static ::ExceptionBase & ExcCellNotActive()
void set_all_manifold_ids(const types::manifold_id) const
void set_all_manifold_ids(const types::manifold_id)
void set_manifold_id(const types::manifold_id) const
@ valid
Iterator points to a valid object.
boost::integer_range< IncrementableType > iota_view
Definition iota_view.h:46
unsigned int manifold_id
Definition types.h:153