Reference documentation for deal.II version GIT relicensing-37-gdcd04c9418 2024-03-02 19:00:02+00:00
\(\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// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 1998 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_tria_accessor_h
16#define dealii_tria_accessor_h
17
18
19#include <deal.II/base/config.h>
20
24#include <deal.II/base/point.h>
25
30
31#include <boost/container/small_vector.hpp>
32
33#include <utility>
34
35
37
38// Forward declarations
39#ifndef DOXYGEN
40template <int dim, int spacedim>
42class Triangulation;
43template <typename Accessor>
44class TriaRawIterator;
45template <typename Accessor>
46class TriaIterator;
47template <typename Accessor>
49
50namespace parallel
51{
52 template <int dim, int spacedim>
54 class TriangulationBase;
55}
56
57template <int dim, int spacedim>
59class DoFHandler;
60template <int dim, int spacedim, bool lda>
61class DoFCellAccessor;
62
63
64template <int dim, int spacedim>
65class Manifold;
66
67template <int dim, int spacedim>
68class Mapping;
69#endif
70
71namespace internal
72{
73 namespace TriangulationImplementation
74 {
75 class TriaObjects;
76 struct Implementation;
77 struct ImplementationMixedMesh;
78 } // namespace TriangulationImplementation
79
80 namespace TriaAccessorImplementation
81 {
82 struct Implementation;
83
89 template <int structdim, int dim>
91 {
92 struct type
93 {
97 type() = default;
98
102 type(const int level)
103 {
105 (void)level; // removes -Wunused-parameter warning in optimized mode
106 }
107
111 operator int() const
112 {
113 return 0;
114 }
115
116 void
118 {
120 }
121
122 void
124 {
126 }
127 };
128 };
129
130
136 template <int dim>
137 struct PresentLevelType<dim, dim>
138 {
139 using type = int;
140 };
141 } // namespace TriaAccessorImplementation
142} // namespace internal
143template <int structdim, int dim, int spacedim>
144class TriaAccessor;
145template <int dim, int spacedim>
146class TriaAccessor<0, dim, spacedim>;
147template <int spacedim>
148class TriaAccessor<0, 1, spacedim>;
149
150// note: the file tria_accessor.templates.h is included at the end of
151// this file. this includes a lot of templates. originally, this was
152// only done in debug mode, but led to cyclic reduction problems and
153// so is now on by default.
154
155
160{
165 "The operation you are attempting can only be performed for "
166 "(cell, face, or edge) iterators that point to valid "
167 "objects. These objects need not necessarily be active, "
168 "i.e., have no children, but they need to be part of a "
169 "triangulation. (The objects pointed to by an iterator "
170 "may -- after coarsening -- also be objects that used "
171 "to be part of a triangulation, but are now no longer "
172 "used. Their memory location may have been retained "
173 "for re-use upon the next mesh refinement, but is "
174 "currently unused.)");
185 "The operation you are attempting can only be performed for "
186 "(cell, face, or edge) iterators that point to 'active' "
187 "objects. 'Active' objects are those that do not have "
188 "children (in the case of cells), or that are part of "
189 "an active cell (in the case of faces or edges). However, "
190 "the object on which you are trying the current "
191 "operation is not 'active' in this sense.");
198 "The operation you are attempting can only be performed for "
199 "(cell, face, or edge) iterators that have children, "
200 "but the object on which you are trying the current "
201 "operation does not have any.");
209 "The operation you are attempting can only be performed for "
210 "(cell, face, or edge) iterators that have a parent object, "
211 "but the object on which you are trying the current "
212 "operation does not have one -- i.e., it is on the "
213 "coarsest level of the triangulation.");
218 int,
219 << "You can only set the child index if the cell does not "
220 << "currently have children registered; or you can clear it. "
221 << "The given index was " << arg1
222 << " (-1 means: clear children).");
226 template <typename AccessorType>
228 AccessorType,
229 << "You tried to dereference an iterator for which this "
230 << "is not possible. More information on this iterator: "
231 << "index=" << arg1.index() << ", state="
232 << (arg1.state() == IteratorState::valid ?
233 "valid" :
234 (arg1.state() == IteratorState::past_the_end ?
235 "past_the_end" :
236 "invalid")));
241 "Iterators can only be compared if they point to the same "
242 "triangulation, or if neither of them are associated "
243 "with a triangulation.");
244 // TODO: Write documentation!
249 // TODO: Write documentation!
269 // TODO: Write documentation!
275 int,
276 << "You can only set the child index of an even numbered child."
277 << "The number of the child given was " << arg1 << '.');
278} // namespace TriaAccessorExceptions
279
280
306template <int structdim, int dim, int spacedim = dim>
308{
309public:
316 static constexpr unsigned int space_dimension = spacedim;
317
323 static constexpr unsigned int dimension = dim;
324
330 static const unsigned int structure_dimension = structdim;
331
341 void
342 operator=(const TriaAccessorBase *) = delete;
343
344protected:
350 using AccessorData = void;
351
356 const int level = -1,
357 const int index = -1,
358 const AccessorData * = nullptr);
359
364
372 void
374
380
391 bool
392 operator<(const TriaAccessorBase &other) const;
393
394protected:
398 bool
400
404 bool
406
420 void
422
430 void
440 objects() const;
441
442public:
448 using LocalData = void *;
449
473 int
474 level() const;
475
502 int
503 index() const;
504
510 state() const;
511
518
522protected:
527 typename ::internal::TriaAccessorImplementation::
528 PresentLevelType<structdim, dim>::type present_level;
529
535
540
541private:
542 template <typename Accessor>
543 friend class TriaRawIterator;
544 template <typename Accessor>
545 friend class TriaIterator;
546 template <typename Accessor>
547 friend class TriaActiveIterator;
548};
549
550
551
572template <int structdim, int dim, int spacedim = dim>
574{
575public:
582 static constexpr unsigned int space_dimension = spacedim;
583
589 static constexpr unsigned int dimension = dim;
590
596 static const unsigned int structure_dimension = structdim;
597
603 using AccessorData = void;
604
612 InvalidAccessor(const void *parent = nullptr,
613 const int level = -1,
614 const int index = -1,
615 const AccessorData *local_data = nullptr);
616
625
630 template <typename OtherAccessor>
631 InvalidAccessor(const OtherAccessor &);
632
636 void
638
642 bool
644 bool
646
650 void
651 operator++() const;
652 void
653 operator--() const;
654
660 state();
661
662
667 static int
668 level();
669
674 static int
675 index();
676
681 bool
682 used() const;
683
688 bool
690
695 manifold_id() const;
696
700 unsigned int
701 user_index() const;
702
706 void
707 set_user_index(const unsigned int p) const;
708
712 void
714
719 vertex(const unsigned int i) const;
720
725 void *
726 line(const unsigned int i) const;
727
732 void *
733 quad(const unsigned int i) const;
734};
735
736
737
755template <int structdim, int dim, int spacedim>
756class TriaAccessor : public TriaAccessorBase<structdim, dim, spacedim>
757{
758public:
764
769 const int level = -1,
770 const int index = -1,
771 const AccessorData *local_data = nullptr);
772
777 TriaAccessor(const TriaAccessor &) = default;
778
782 TriaAccessor(TriaAccessor &&) = default; // NOLINT
783
796 template <int structdim2, int dim2, int spacedim2>
798
803 template <int structdim2, int dim2, int spacedim2>
805
816 operator=(const TriaAccessor &) = delete;
817
822 operator=(TriaAccessor &&) = default; // NOLINT
823
827 ~TriaAccessor() = default;
828
835 bool
836 used() const;
837
850 vertex_iterator(const unsigned int i) const;
851
867 unsigned int
868 vertex_index(const unsigned int i) const;
869
908 vertex(const unsigned int i) const;
909
913 typename ::internal::TriangulationImplementation::
914 Iterators<dim, spacedim>::line_iterator
915 line(const unsigned int i) const;
916
923 unsigned int
924 line_index(const unsigned int i) const;
925
929 typename ::internal::TriangulationImplementation::
930 Iterators<dim, spacedim>::quad_iterator
931 quad(const unsigned int i) const;
932
939 unsigned int
940 quad_index(const unsigned int i) const;
958 unsigned char
959 combined_face_orientation(const unsigned int face) const;
960
972 bool
973 face_orientation(const unsigned int face) const;
974
984 bool
985 face_flip(const unsigned int face) const;
986
996 bool
997 face_rotation(const unsigned int face) const;
998
1011 bool
1012 line_orientation(const unsigned int line) const;
1027 bool
1029
1034 unsigned int
1035 n_children() const;
1036
1050 unsigned int
1052
1066 unsigned int
1068
1073 child(const unsigned int i) const;
1074
1079 unsigned int
1082
1092 isotropic_child(const unsigned int i) const;
1093
1099
1105 int
1106 child_index(const unsigned int i) const;
1107
1113 int
1114 isotropic_child_index(const unsigned int i) const;
1138
1168 void
1170
1201 void
1203
1211 bool
1213
1225
1248
1266 void
1268
1282 void
1284
1301 bool
1303
1309 void
1311
1317 void
1319
1325 void
1327
1333 void
1335
1341 void
1343
1355 void
1356 set_user_pointer(void *p) const;
1357
1363 void
1365
1381 void *
1383
1405 void
1407
1414 void
1416
1426 void
1427 set_user_index(const unsigned int p) const;
1428
1434 void
1436
1448 unsigned int
1449 user_index() const;
1450
1468 void
1469 recursively_set_user_index(const unsigned int p) const;
1470
1479 void
1515 double
1516 diameter() const;
1517
1544 std::pair<Point<spacedim>, double>
1546
1557
1567 double
1568 extent_in_direction(const unsigned int axis) const;
1569
1573 double
1575
1590 intermediate_point(const Point<structdim> &coordinates) const;
1591
1616
1652 center(const bool respect_manifold = false,
1653 const bool interpolate_from_surrounding = false) const;
1654
1673 barycenter() const;
1674
1705 double
1706 measure() const;
1707
1722 bool
1725
1731
1735 unsigned int
1736 n_vertices() const;
1737
1741 unsigned int
1742 n_lines() const;
1743
1753 unsigned int
1754 n_faces() const;
1755
1762
1769
1778
1783private:
1788 void
1790
1798 void
1800 const std::initializer_list<int> &new_indices) const;
1801
1805 void
1807 const std::initializer_list<unsigned int> &new_indices) const;
1808
1816 void
1817 set_line_orientation(const unsigned int line, const bool orientation) const;
1818
1828 void
1829 set_combined_face_orientation(const unsigned int face,
1830 const unsigned char combined_orientation) const;
1831
1835 void
1837
1841 void
1843
1852 void
1854
1862 void
1864
1871 void
1872 set_children(const unsigned int i, const int index) const;
1873
1878 void
1880
1881private:
1882 friend class Triangulation<dim, spacedim>;
1883
1884 friend struct ::internal::TriangulationImplementation::Implementation;
1885 friend struct ::internal::TriangulationImplementation::
1886 ImplementationMixedMesh;
1887 friend struct ::internal::TriaAccessorImplementation::Implementation;
1888};
1889
1890
1891
1910template <int dim, int spacedim>
1911class TriaAccessor<0, dim, spacedim>
1912{
1913public:
1919 static constexpr unsigned int space_dimension = spacedim;
1920
1926 static constexpr unsigned int dimension = dim;
1927
1933 static const unsigned int structure_dimension = 0;
1934
1938 using AccessorData = void;
1939
1945 const unsigned int vertex_index);
1946
1953 const int level = 0,
1954 const int index = 0,
1955 const AccessorData * = nullptr);
1956
1960 template <int structdim2, int dim2, int spacedim2>
1962
1966 template <int structdim2, int dim2, int spacedim2>
1968
1973 state() const;
1974
1979 static int
1981
1986 int
1987 index() const;
1988
1995
2005 void
2007
2011 void
2016 bool
2017 operator==(const TriaAccessor &) const;
2018
2022 bool
2023 operator!=(const TriaAccessor &) const;
2024
2052 unsigned int
2053 vertex_index(const unsigned int i = 0) const;
2054
2061 vertex(const unsigned int i = 0) const;
2062
2067 typename ::internal::TriangulationImplementation::
2068 Iterators<dim, spacedim>::line_iterator static line(const unsigned int);
2069
2073 static unsigned int
2074 line_index(const unsigned int i);
2075
2079 static typename ::internal::TriangulationImplementation::
2080 Iterators<dim, spacedim>::quad_iterator
2081 quad(const unsigned int i);
2082
2086 static unsigned int
2087 quad_index(const unsigned int i);
2088
2104 double
2105 diameter() const;
2106
2114 double
2115 extent_in_direction(const unsigned int axis) const;
2116
2125 center(const bool respect_manifold = false,
2126 const bool interpolate_from_surrounding = false) const;
2127
2136 double
2137 measure() const;
2154 static unsigned char
2155 combined_face_orientation(const unsigned int face);
2156
2160 static bool
2161 face_orientation(const unsigned int face);
2162
2166 static bool
2167 face_flip(const unsigned int face);
2168
2172 static bool
2173 face_rotation(const unsigned int face);
2174
2178 static bool
2179 line_orientation(const unsigned int line);
2180
2195 static bool
2197
2202 static unsigned int
2204
2209 static unsigned int
2211
2215 static unsigned int
2217
2221 static unsigned int
2223
2228 child(const unsigned int);
2229
2234 isotropic_child(const unsigned int);
2235
2239 static RefinementCase<0>
2241
2245 static int
2246 child_index(const unsigned int i);
2247
2251 static int
2252 isotropic_child_index(const unsigned int i);
2260 bool
2261 used() const;
2262
2263protected:
2271 void
2273
2282 bool
2283 operator<(const TriaAccessor &other) const;
2284
2289
2294
2295private:
2296 template <typename Accessor>
2297 friend class TriaRawIterator;
2298 template <typename Accessor>
2299 friend class TriaIterator;
2300 template <typename Accessor>
2302};
2303
2304
2305
2322template <int spacedim>
2323class TriaAccessor<0, 1, spacedim>
2324{
2325public:
2331 static constexpr unsigned int space_dimension = spacedim;
2332
2338 static constexpr unsigned int dimension = 1;
2339
2345 static const unsigned int structure_dimension = 0;
2346
2350 using AccessorData = void;
2351
2357 {
2369 right_vertex
2371
2384 const VertexKind vertex_kind,
2385 const unsigned int vertex_index);
2386
2393 const int = 0,
2394 const int = 0,
2395 const AccessorData * = nullptr);
2396
2400 template <int structdim2, int dim2, int spacedim2>
2402
2406 template <int structdim2, int dim2, int spacedim2>
2408
2413 void
2415
2421 void
2423
2431
2436 static int
2438
2443 int
2444 index() const;
2445
2452
2463 void
2464 operator++() const;
2465
2470 void
2471 operator--() const;
2475 bool
2476 operator==(const TriaAccessor &) const;
2477
2481 bool
2482 operator!=(const TriaAccessor &) const;
2483
2492 bool
2493 operator<(const TriaAccessor &other) const;
2494
2521 unsigned int
2522 vertex_index(const unsigned int i = 0) const;
2523
2530 vertex(const unsigned int i = 0) const;
2531
2537 center() const;
2538
2543 typename ::internal::TriangulationImplementation::
2544 Iterators<1, spacedim>::line_iterator static line(const unsigned int);
2545
2552 static unsigned int
2553 line_index(const unsigned int i);
2554
2558 static typename ::internal::TriangulationImplementation::
2559 Iterators<1, spacedim>::quad_iterator
2560 quad(const unsigned int i);
2561
2568 static unsigned int
2569 quad_index(const unsigned int i);
2570
2580 bool
2582
2599
2603 const Manifold<1, spacedim> &
2605
2614
2615
2627 bool
2628 user_flag_set() const;
2629
2635 void
2636 set_user_flag() const;
2637
2643 void
2644 clear_user_flag() const;
2645
2651 void
2653
2659 void
2661
2667 void
2668 clear_user_data() const;
2669
2681 void
2682 set_user_pointer(void *p) const;
2683
2689 void
2690 clear_user_pointer() const;
2691
2707 void *
2708 user_pointer() const;
2709
2731 void
2732 recursively_set_user_pointer(void *p) const;
2733
2740 void
2742
2752 void
2753 set_user_index(const unsigned int p) const;
2754
2760 void
2761 clear_user_index() const;
2762
2774 unsigned int
2775 user_index() const;
2776
2794 void
2795 recursively_set_user_index(const unsigned int p) const;
2796
2805 void
2823 static unsigned char
2824 combined_face_orientation(const unsigned int face);
2825
2829 static bool
2830 face_orientation(const unsigned int face);
2831
2835 static bool
2836 face_flip(const unsigned int face);
2837
2841 static bool
2842 face_rotation(const unsigned int face);
2843
2847 static bool
2848 line_orientation(const unsigned int line);
2849
2864 static bool
2866
2871 static unsigned int
2873
2878 static unsigned int
2880
2884 static unsigned int
2886
2890 static unsigned int
2892
2897 child(const unsigned int);
2898
2903 isotropic_child(const unsigned int);
2904
2908 static RefinementCase<0>
2910
2914 static int
2915 child_index(const unsigned int i);
2916
2920 static int
2921 isotropic_child_index(const unsigned int i);
2954 void
2956
2963 void
2965
2977 void
2979
2991 void
3000 bool
3001 used() const;
3002
3008
3012 unsigned int
3013 n_vertices() const;
3014
3018 unsigned int
3019 n_lines() const;
3020
3027
3034
3035protected:
3040
3046
3051};
3052
3053
3054
3070template <int dim, int spacedim = dim>
3071class CellAccessor : public TriaAccessor<dim, dim, spacedim>
3072{
3073public:
3078
3083
3095 const int level = -1,
3096 const int index = -1,
3097 const AccessorData *local_data = nullptr);
3098
3103
3116 template <int structdim2, int dim2, int spacedim2>
3118
3123 template <int structdim2, int dim2, int spacedim2>
3125
3130
3134 // NOLINTNEXTLINE OSX does not compile with noexcept
3136
3140 ~CellAccessor() = default;
3141
3153
3157 // NOLINTNEXTLINE OSX does not compile with noexcept
3160
3184 as_dof_handler_iterator(const DoFHandler<dim, spacedim> &dof_handler) const;
3185
3196 const DoFHandler<dim, spacedim> &dof_handler) const;
3197
3198
3215 child(const unsigned int i) const;
3216
3220 boost::container::small_vector<TriaIterator<CellAccessor<dim, spacedim>>,
3223
3227 TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>
3228 face(const unsigned int i) const;
3229
3234 unsigned int
3237
3241 boost::container::small_vector<
3242 TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>,
3245
3255 unsigned int
3256 face_index(const unsigned int i) const;
3257
3306 neighbor_child_on_subface(const unsigned int face_no,
3307 const unsigned int subface_no) const;
3308
3334 neighbor(const unsigned int face_no) const;
3335
3343 int
3344 neighbor_index(const unsigned int face_no) const;
3345
3353 int
3354 neighbor_level(const unsigned int face_no) const;
3355
3367 unsigned int
3368 neighbor_of_neighbor(const unsigned int face_no) const;
3369
3380 bool
3381 neighbor_is_coarser(const unsigned int face_no) const;
3382
3397 std::pair<unsigned int, unsigned int>
3398 neighbor_of_coarser_neighbor(const unsigned int neighbor) const;
3399
3406 unsigned int
3407 neighbor_face_no(const unsigned int neighbor) const;
3408
3412 static bool
3414
3428 bool
3429 has_periodic_neighbor(const unsigned int i) const;
3430
3448 periodic_neighbor(const unsigned int i) const;
3449
3458 neighbor_or_periodic_neighbor(const unsigned int i) const;
3459
3475 periodic_neighbor_child_on_subface(const unsigned int face_no,
3476 const unsigned int subface_no) const;
3477
3488 std::pair<unsigned int, unsigned int>
3489 periodic_neighbor_of_coarser_periodic_neighbor(const unsigned face_no) const;
3490
3496 int
3497 periodic_neighbor_index(const unsigned int i) const;
3498
3504 int
3505 periodic_neighbor_level(const unsigned int i) const;
3506
3521 unsigned int
3522 periodic_neighbor_of_periodic_neighbor(const unsigned int i) const;
3523
3529 unsigned int
3530 periodic_neighbor_face_no(const unsigned int i) const;
3531
3538 bool
3539 periodic_neighbor_is_coarser(const unsigned int i) const;
3540
3557 bool
3558 at_boundary(const unsigned int i) const;
3559
3568 bool
3569 at_boundary() const;
3570
3578 bool
3579 has_boundary_lines() const;
3607
3625 void
3628
3632 void
3634
3642 bool
3644 const unsigned int face_no,
3645 const RefinementCase<dim - 1> &face_refinement_case =
3647
3653 bool
3654 flag_for_line_refinement(const unsigned int line_no) const;
3655
3665 subface_case(const unsigned int face_no) const;
3666
3670 bool
3672
3677 void
3679
3683 void
3708 material_id() const;
3709
3721 void
3722 set_material_id(const types::material_id new_material_id) const;
3723
3732 void
3733 recursively_set_material_id(const types::material_id new_material_id) const;
3761
3777 void
3778 set_subdomain_id(const types::subdomain_id new_subdomain_id) const;
3779
3790
3795 void
3797 const types::subdomain_id new_level_subdomain_id) const;
3798
3799
3815 void
3817 const types::subdomain_id new_subdomain_id) const;
3841
3851
3868 bool
3869 direction_flag() const;
3870
3896 unsigned int
3898
3906 int
3907 parent_index() const;
3908
3915 parent() const;
3916
3936 bool
3937 is_active() const;
3938
3958 bool
3960
3965 bool
3967
4001 bool
4002 is_ghost() const;
4003
4009 bool
4011
4038 bool
4040
4047 bool
4049
4063 bool
4065
4074 void
4075 set_neighbor(const unsigned int i,
4076 const TriaIterator<CellAccessor<dim, spacedim>> &pointer) const;
4077
4091 CellId
4092 id() const;
4093
4094 using TriaAccessor<dim, dim, spacedim>::diameter;
4095
4099 double
4100 diameter(const Mapping<dim, spacedim> &mapping) const;
4101
4119
4120protected:
4136 unsigned int
4137 neighbor_of_neighbor_internal(const unsigned int neighbor) const;
4138
4144 template <int dim_, int spacedim_>
4145 bool
4146 point_inside_codim(const Point<spacedim_> &p) const;
4147
4148
4149
4150private:
4155 void
4156 set_active_cell_index(const unsigned int active_cell_index) const;
4157
4161 void
4163
4167 void
4169
4173 void
4174 set_parent(const unsigned int parent_index);
4175
4186 void
4187 set_direction_flag(const bool new_direction_flag) const;
4188
4189 friend class Triangulation<dim, spacedim>;
4190
4191 friend class parallel::TriangulationBase<dim, spacedim>;
4192
4193 friend struct ::internal::TriangulationImplementation::Implementation;
4194 friend struct ::internal::TriangulationImplementation::
4195 ImplementationMixedMesh;
4196};
4197
4198
4199
4200/* ----- declaration of explicit specializations and general templates ----- */
4201
4202
4203template <int structdim, int dim, int spacedim>
4204template <typename OtherAccessor>
4206 const OtherAccessor &)
4207{
4208 Assert(false,
4209 ExcMessage("You are attempting an illegal conversion between "
4210 "iterator/accessor types. The constructor you call "
4211 "only exists to make certain template constructs "
4212 "easier to write as dimension independent code but "
4213 "the conversion is not valid in the current context."));
4214}
4215
4216
4217
4218template <int structdim, int dim, int spacedim>
4219template <int structdim2, int dim2, int spacedim2>
4222{
4223 Assert(false,
4224 ExcMessage("You are attempting an illegal conversion between "
4225 "iterator/accessor types. The constructor you call "
4226 "only exists to make certain template constructs "
4227 "easier to write as dimension independent code but "
4228 "the conversion is not valid in the current context."));
4229}
4230
4231
4232
4233template <int dim, int spacedim>
4234template <int structdim2, int dim2, int spacedim2>
4237{
4238 Assert(false,
4239 ExcMessage("You are attempting an illegal conversion between "
4240 "iterator/accessor types. The constructor you call "
4241 "only exists to make certain template constructs "
4242 "easier to write as dimension independent code but "
4243 "the conversion is not valid in the current context."));
4244}
4245
4246
4247
4248template <int structdim, int dim, int spacedim>
4249template <int structdim2, int dim2, int spacedim2>
4252{
4253 Assert(false,
4254 ExcMessage("You are attempting an illegal conversion between "
4255 "iterator/accessor types. The constructor you call "
4256 "only exists to make certain template constructs "
4257 "easier to write as dimension independent code but "
4258 "the conversion is not valid in the current context."));
4259}
4260
4261
4262
4263template <int dim, int spacedim>
4264template <int structdim2, int dim2, int spacedim2>
4267{
4268 Assert(false,
4269 ExcMessage("You are attempting an illegal conversion between "
4270 "iterator/accessor types. The constructor you call "
4271 "only exists to make certain template constructs "
4272 "easier to write as dimension independent code but "
4273 "the conversion is not valid in the current context."));
4274}
4275
4276
4277#ifndef DOXYGEN
4278
4279template <>
4280bool
4282template <>
4283bool
4285template <>
4286bool
4288template <>
4289bool
4291template <>
4292bool
4294template <>
4295bool
4297// -------------------------------------------------------------------
4298
4299template <>
4300void
4302
4303#endif // DOXYGEN
4304
4306
4307// include more templates in debug and optimized mode
4308#include <deal.II/grid/tria_accessor.templates.h>
4309
4310#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
TriaIterator< DoFCellAccessor< dim, spacedim, true > > as_dof_handler_level_iterator(const DoFHandler< dim, spacedim > &dof_handler) 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:316
Definition point.h:111
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 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 > &)
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 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()
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 > &)
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 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 unsigned int max_refinement_depth()
static unsigned int quad_index(const unsigned int i)
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
std_cxx20::ranges::iota_view< unsigned int, unsigned int > vertex_indices() const
void clear_children() 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
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
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 set_bounding_object_indices(const std::initializer_list< unsigned int > &new_indices) 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 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
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
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_CXX20_REQUIRES(condition)
Definition config.h:177
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
Point< 3 > center
unsigned int level
Definition grid_out.cc:4616
static ::ExceptionBase & ExcCantCompareIterators()
#define DeclException0(Exception0)
Definition exceptions.h:471
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:494
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcFacesHaveNoLevel()
static ::ExceptionBase & ExcDereferenceInvalidObject(AccessorType arg1)
static ::ExceptionBase & ExcCellFlaggedForCoarsening()
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:516
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
bool line_orientation(const unsigned int line) const
static bool face_orientation(const unsigned int face)
Always return false.
unsigned char combined_face_orientation(const unsigned int face) const
bool face_rotation(const unsigned int face) const
bool face_orientation(const unsigned int face) const
static bool face_flip(const unsigned int face)
Always return false.
static bool line_orientation(const unsigned int line)
Always return false.
static unsigned char combined_face_orientation(const unsigned int face)
Always return 0.
static bool face_flip(const unsigned int face)
Always return false.
static bool line_orientation(const unsigned int line)
Always return false.
static bool face_rotation(const unsigned int face)
Always return false.
static bool face_rotation(const unsigned int face)
Always return false.
static bool face_orientation(const unsigned int face)
Always return false.
static unsigned char combined_face_orientation(const unsigned int face)
Always return 0.
bool face_flip(const unsigned int face) const
void set_combined_face_orientation(const unsigned int face, const unsigned char combined_orientation) const
#define DEAL_II_ASSERT_UNREACHABLE()
@ valid
Iterator points to a valid object.
boost::integer_range< IncrementableType > iota_view
Definition iota_view.h:45