Reference documentation for deal.II version GIT relicensing-397-g31a1263477 2024-04-16 03:30: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
1826 void
1827 set_combined_face_orientation(const unsigned int face,
1828 const unsigned char combined_orientation) const;
1829
1833 void
1835
1839 void
1841
1850 void
1852
1860 void
1862
1869 void
1870 set_children(const unsigned int i, const int index) const;
1871
1876 void
1878
1879private:
1880 friend class Triangulation<dim, spacedim>;
1881
1882 friend struct ::internal::TriangulationImplementation::Implementation;
1883 friend struct ::internal::TriangulationImplementation::
1884 ImplementationMixedMesh;
1885 friend struct ::internal::TriaAccessorImplementation::Implementation;
1886};
1887
1888
1889
1908template <int dim, int spacedim>
1909class TriaAccessor<0, dim, spacedim>
1910{
1911public:
1917 static constexpr unsigned int space_dimension = spacedim;
1918
1924 static constexpr unsigned int dimension = dim;
1925
1931 static const unsigned int structure_dimension = 0;
1932
1936 using AccessorData = void;
1937
1943 const unsigned int vertex_index);
1944
1951 const int level = 0,
1952 const int index = 0,
1953 const AccessorData * = nullptr);
1954
1958 template <int structdim2, int dim2, int spacedim2>
1960
1964 template <int structdim2, int dim2, int spacedim2>
1966
1971 state() const;
1972
1977 static int
1979
1984 int
1985 index() const;
1986
1993
2003 void
2005
2009 void
2014 bool
2015 operator==(const TriaAccessor &) const;
2016
2020 bool
2021 operator!=(const TriaAccessor &) const;
2022
2050 unsigned int
2051 vertex_index(const unsigned int i = 0) const;
2052
2059 vertex(const unsigned int i = 0) const;
2060
2065 typename ::internal::TriangulationImplementation::
2066 Iterators<dim, spacedim>::line_iterator static line(const unsigned int);
2067
2071 static unsigned int
2072 line_index(const unsigned int i);
2073
2077 static typename ::internal::TriangulationImplementation::
2078 Iterators<dim, spacedim>::quad_iterator
2079 quad(const unsigned int i);
2080
2084 static unsigned int
2085 quad_index(const unsigned int i);
2086
2102 double
2103 diameter() const;
2104
2112 double
2113 extent_in_direction(const unsigned int axis) const;
2114
2123 center(const bool respect_manifold = false,
2124 const bool interpolate_from_surrounding = false) const;
2125
2134 double
2135 measure() const;
2152 static unsigned char
2153 combined_face_orientation(const unsigned int face);
2154
2158 static bool
2159 face_orientation(const unsigned int face);
2160
2164 static bool
2165 face_flip(const unsigned int face);
2166
2170 static bool
2171 face_rotation(const unsigned int face);
2172
2176 static bool
2177 line_orientation(const unsigned int line);
2178
2193 static bool
2195
2200 static unsigned int
2202
2207 static unsigned int
2209
2213 static unsigned int
2215
2219 static unsigned int
2221
2226 child(const unsigned int);
2227
2232 isotropic_child(const unsigned int);
2233
2237 static RefinementCase<0>
2239
2243 static int
2244 child_index(const unsigned int i);
2245
2249 static int
2250 isotropic_child_index(const unsigned int i);
2258 bool
2259 used() const;
2260
2261protected:
2269 void
2271
2280 bool
2281 operator<(const TriaAccessor &other) const;
2282
2287
2292
2293private:
2294 template <typename Accessor>
2295 friend class TriaRawIterator;
2296 template <typename Accessor>
2297 friend class TriaIterator;
2298 template <typename Accessor>
2300};
2301
2302
2303
2320template <int spacedim>
2321class TriaAccessor<0, 1, spacedim>
2322{
2323public:
2329 static constexpr unsigned int space_dimension = spacedim;
2330
2336 static constexpr unsigned int dimension = 1;
2337
2343 static const unsigned int structure_dimension = 0;
2344
2348 using AccessorData = void;
2349
2355 {
2367 right_vertex
2369
2382 const VertexKind vertex_kind,
2383 const unsigned int vertex_index);
2384
2391 const int = 0,
2392 const int = 0,
2393 const AccessorData * = nullptr);
2394
2398 template <int structdim2, int dim2, int spacedim2>
2400
2404 template <int structdim2, int dim2, int spacedim2>
2406
2411 void
2413
2419 void
2421
2429
2434 static int
2436
2441 int
2442 index() const;
2443
2450
2461 void
2462 operator++() const;
2463
2468 void
2469 operator--() const;
2473 bool
2474 operator==(const TriaAccessor &) const;
2475
2479 bool
2480 operator!=(const TriaAccessor &) const;
2481
2490 bool
2491 operator<(const TriaAccessor &other) const;
2492
2519 unsigned int
2520 vertex_index(const unsigned int i = 0) const;
2521
2528 vertex(const unsigned int i = 0) const;
2529
2535 center() const;
2536
2541 typename ::internal::TriangulationImplementation::
2542 Iterators<1, spacedim>::line_iterator static line(const unsigned int);
2543
2550 static unsigned int
2551 line_index(const unsigned int i);
2552
2556 static typename ::internal::TriangulationImplementation::
2557 Iterators<1, spacedim>::quad_iterator
2558 quad(const unsigned int i);
2559
2566 static unsigned int
2567 quad_index(const unsigned int i);
2568
2578 bool
2580
2597
2601 const Manifold<1, spacedim> &
2603
2612
2613
2625 bool
2626 user_flag_set() const;
2627
2633 void
2634 set_user_flag() const;
2635
2641 void
2642 clear_user_flag() const;
2643
2649 void
2651
2657 void
2659
2665 void
2666 clear_user_data() const;
2667
2679 void
2680 set_user_pointer(void *p) const;
2681
2687 void
2688 clear_user_pointer() const;
2689
2705 void *
2706 user_pointer() const;
2707
2729 void
2730 recursively_set_user_pointer(void *p) const;
2731
2738 void
2740
2750 void
2751 set_user_index(const unsigned int p) const;
2752
2758 void
2759 clear_user_index() const;
2760
2772 unsigned int
2773 user_index() const;
2774
2792 void
2793 recursively_set_user_index(const unsigned int p) const;
2794
2803 void
2821 static unsigned char
2822 combined_face_orientation(const unsigned int face);
2823
2827 static bool
2828 face_orientation(const unsigned int face);
2829
2833 static bool
2834 face_flip(const unsigned int face);
2835
2839 static bool
2840 face_rotation(const unsigned int face);
2841
2845 static bool
2846 line_orientation(const unsigned int line);
2847
2862 static bool
2864
2869 static unsigned int
2871
2876 static unsigned int
2878
2882 static unsigned int
2884
2888 static unsigned int
2890
2895 child(const unsigned int);
2896
2901 isotropic_child(const unsigned int);
2902
2906 static RefinementCase<0>
2908
2912 static int
2913 child_index(const unsigned int i);
2914
2918 static int
2919 isotropic_child_index(const unsigned int i);
2952 void
2954
2961 void
2963
2975 void
2977
2989 void
2998 bool
2999 used() const;
3000
3006
3010 unsigned int
3011 n_vertices() const;
3012
3016 unsigned int
3017 n_lines() const;
3018
3025
3032
3033protected:
3038
3044
3049};
3050
3051
3052
3068template <int dim, int spacedim = dim>
3069class CellAccessor : public TriaAccessor<dim, dim, spacedim>
3070{
3071public:
3076
3081
3093 const int level = -1,
3094 const int index = -1,
3095 const AccessorData *local_data = nullptr);
3096
3101
3114 template <int structdim2, int dim2, int spacedim2>
3116
3121 template <int structdim2, int dim2, int spacedim2>
3123
3128
3132 // NOLINTNEXTLINE OSX does not compile with noexcept
3134
3138 ~CellAccessor() = default;
3139
3151
3155 // NOLINTNEXTLINE OSX does not compile with noexcept
3158
3182 as_dof_handler_iterator(const DoFHandler<dim, spacedim> &dof_handler) const;
3183
3194 const DoFHandler<dim, spacedim> &dof_handler) const;
3195
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
3866 bool
3867 direction_flag() const;
3868
3894 unsigned int
3896
3904 int
3905 parent_index() const;
3906
3913 parent() const;
3914
3934 bool
3935 is_active() const;
3936
3956 bool
3958
3963 bool
3965
3999 bool
4000 is_ghost() const;
4001
4007 bool
4009
4036 bool
4038
4045 bool
4047
4061 bool
4063
4072 void
4073 set_neighbor(const unsigned int i,
4074 const TriaIterator<CellAccessor<dim, spacedim>> &pointer) const;
4075
4089 CellId
4090 id() const;
4091
4092 using TriaAccessor<dim, dim, spacedim>::diameter;
4093
4097 double
4098 diameter(const Mapping<dim, spacedim> &mapping) const;
4099
4117
4118protected:
4134 unsigned int
4135 neighbor_of_neighbor_internal(const unsigned int neighbor) const;
4136
4142 template <int dim_, int spacedim_>
4143 bool
4144 point_inside_codim(const Point<spacedim_> &p) const;
4145
4146
4147
4148private:
4153 void
4154 set_active_cell_index(const unsigned int active_cell_index) const;
4155
4159 void
4161
4165 void
4167
4171 void
4172 set_parent(const unsigned int parent_index);
4173
4184 void
4185 set_direction_flag(const bool new_direction_flag) const;
4186
4187 friend class Triangulation<dim, spacedim>;
4188
4189 friend class parallel::TriangulationBase<dim, spacedim>;
4190
4191 friend struct ::internal::TriangulationImplementation::Implementation;
4192 friend struct ::internal::TriangulationImplementation::
4193 ImplementationMixedMesh;
4194};
4195
4196
4197
4198/* ----- declaration of explicit specializations and general templates ----- */
4199
4200
4201template <int structdim, int dim, int spacedim>
4202template <typename OtherAccessor>
4204 const OtherAccessor &)
4205{
4206 Assert(false,
4207 ExcMessage("You are attempting an illegal conversion between "
4208 "iterator/accessor types. The constructor you call "
4209 "only exists to make certain template constructs "
4210 "easier to write as dimension independent code but "
4211 "the conversion is not valid in the current context."));
4212}
4213
4214
4215
4216template <int structdim, int dim, int spacedim>
4217template <int structdim2, int dim2, int spacedim2>
4220{
4221 Assert(false,
4222 ExcMessage("You are attempting an illegal conversion between "
4223 "iterator/accessor types. The constructor you call "
4224 "only exists to make certain template constructs "
4225 "easier to write as dimension independent code but "
4226 "the conversion is not valid in the current context."));
4227}
4228
4229
4230
4231template <int dim, int spacedim>
4232template <int structdim2, int dim2, int spacedim2>
4235{
4236 Assert(false,
4237 ExcMessage("You are attempting an illegal conversion between "
4238 "iterator/accessor types. The constructor you call "
4239 "only exists to make certain template constructs "
4240 "easier to write as dimension independent code but "
4241 "the conversion is not valid in the current context."));
4242}
4243
4244
4245
4246template <int structdim, int dim, int spacedim>
4247template <int structdim2, int dim2, int spacedim2>
4250{
4251 Assert(false,
4252 ExcMessage("You are attempting an illegal conversion between "
4253 "iterator/accessor types. The constructor you call "
4254 "only exists to make certain template constructs "
4255 "easier to write as dimension independent code but "
4256 "the conversion is not valid in the current context."));
4257}
4258
4259
4260
4261template <int dim, int spacedim>
4262template <int structdim2, int dim2, int spacedim2>
4265{
4266 Assert(false,
4267 ExcMessage("You are attempting an illegal conversion between "
4268 "iterator/accessor types. The constructor you call "
4269 "only exists to make certain template constructs "
4270 "easier to write as dimension independent code but "
4271 "the conversion is not valid in the current context."));
4272}
4273
4274
4275#ifndef DOXYGEN
4276
4277template <>
4278bool
4280template <>
4281bool
4283template <>
4284bool
4286template <>
4287bool
4289template <>
4290bool
4292template <>
4293bool
4295// -------------------------------------------------------------------
4296
4297template <>
4298void
4300
4301#endif // DOXYGEN
4302
4304
4305// include more templates in debug and optimized mode
4306#include <deal.II/grid/tria_accessor.templates.h>
4307
4308#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:318
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:4626
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