15#ifndef dealii_dof_accessor_h
16#define dealii_dof_accessor_h
33#include <boost/container/small_vector.hpp>
45template <
typename number>
47template <
typename number>
49template <
typename number>
52template <
typename Accessor>
60 namespace DoFCellAccessorImplementation
62 struct Implementation;
65 namespace DoFHandlerImplementation
67 struct Implementation;
70 struct Implementation;
76 namespace DoFHandlerImplementation
78 struct Implementation;
87 namespace DoFAccessorImplementation
105 template <
int structdim,
int dim,
int spacedim>
121 template <
int dim,
int spacedim>
87 namespace DoFAccessorImplementation {
…}
212template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
233 using BaseClass = typename ::internal::DoFAccessorImplementation::
234 Inheritance<structdim, dimension, space_dimension>::BaseClass;
303 template <
int structdim2,
int dim2,
int spacedim2>
310 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
317 template <
bool level_dof_access2>
354 template <
bool level_dof_access2>
390 typename ::internal::DoFHandlerImplementation::
391 Iterators<dim, spacedim, level_dof_access>::line_iterator
392 line(
const unsigned int i)
const;
399 typename ::internal::DoFHandlerImplementation::
400 Iterators<dim, spacedim, level_dof_access>::quad_iterator
401 quad(
const unsigned int i)
const;
450 std::vector<types::global_dof_index> &dof_indices,
462 std::vector<types::global_dof_index> &dof_indices,
471 const std::vector<types::global_dof_index> &dof_indices,
497 const unsigned int vertex,
498 const unsigned int i,
509 const unsigned int vertex,
510 const unsigned int i,
590 std::set<types::fe_index>
623 "This accessor object has not been "
624 "associated with any DoFHandler object.");
678 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
686 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
718 const unsigned int i,
724 const unsigned int i,
730 const unsigned int vertex,
731 const unsigned int i,
739 template <
int,
int,
int,
bool>
747 friend struct ::internal::DoFHandlerImplementation::Policy::
749 friend struct ::internal::DoFHandlerImplementation::Implementation;
750 friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
751 friend struct ::internal::DoFCellAccessorImplementation::Implementation;
752 friend struct ::internal::DoFAccessorImplementation::Implementation;
764template <
int spacedim,
bool level_dof_access>
823 const unsigned int vertex_index,
849 template <
int structdim2,
int dim2,
int spacedim2>
856 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
908 template <
bool level_dof_access2>
932 child(const
unsigned int c) const;
940 typename ::
internal::DoFHandlerImplementation::
941 Iterators<1, spacedim, level_dof_access>::line_iterator
942 line(const
unsigned int i) const;
950 typename ::
internal::DoFHandlerImplementation::
951 Iterators<1, spacedim, level_dof_access>::quad_iterator
952 quad(const
unsigned int i) const;
999 std::vector<
types::global_dof_index> &dof_indices,
1000 const
types::fe_index fe_index =
numbers::invalid_fe_index) const;
1011 std::vector<
types::global_dof_index> &dof_indices,
1012 const
types::fe_index fe_index =
numbers::invalid_fe_index) const;
1032 types::global_dof_index
1034 const
unsigned int vertex,
1035 const
unsigned int i,
1036 const
types::fe_index fe_index =
numbers::invalid_fe_index) const;
1055 types::global_dof_index
1057 const
types::fe_index fe_index =
numbers::invalid_fe_index) const;
1119 "This accessor
object has not been "
1163 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
1166 const
DoFAccessor<structdim2, dim2, spacedim2, level_dof_access2> &) const;
1171 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
1174 const
DoFAccessor<structdim2, dim2, spacedim2, level_dof_access2> &) const;
1202 const
unsigned int i,
1203 const
types::global_dof_index index,
1204 const
types::fe_index fe_index =
numbers::invalid_fe_index) const;
1216 friend struct ::
internal::DoFHandlerImplementation::Policy::
1248template <
int structdim,
int dim,
int spacedim = dim>
1266 const int level = -1,
1267 const int index = -1,
1283 template <
typename OtherAccessor>
1303 const unsigned int i,
1324template <
int dimension_,
int space_dimension_,
bool level_dof_access>
1334 static const unsigned int dim = dimension_;
1339 static const unsigned int spacedim = space_dimension_;
1395 template <
int structdim2,
int dim2,
int spacedim2>
1402 template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
1476 neighbor(
const unsigned int i)
const;
1484 periodic_neighbor(
const unsigned int i)
const;
1492 neighbor_or_periodic_neighbor(
const unsigned int i)
const;
1500 child(
const unsigned int i)
const;
1505 boost::container::small_vector<
1509 child_iterators()
const;
1518 face(
const unsigned int i)
const;
1525 face_iterators()
const;
1534 neighbor_child_on_subface(
const unsigned int face_no,
1535 const unsigned int subface_no)
const;
1544 periodic_neighbor_child_on_subface(
const unsigned int face_no,
1545 const unsigned int subface_no)
const;
1575 template <
class InputVector,
typename number>
1577 get_dof_values(
const InputVector &values,
Vector<number> &local_values)
const;
1596 template <
typename Number,
typename ForwardIterator>
1599 ForwardIterator local_values_begin,
1600 ForwardIterator local_values_end)
const;
1622 template <
class InputVector,
typename ForwardIterator>
1626 const InputVector &values,
1627 ForwardIterator local_values_begin,
1628 ForwardIterator local_values_end)
const;
1654 template <
class OutputVector,
typename number>
1657 OutputVector &values)
const;
1690 template <
typename Number>
1692 get_interpolated_dof_values(
1758 template <
class OutputVector,
typename number>
1760 set_dof_values_by_interpolation(
1762 OutputVector &values,
1764 const bool perform_check =
false)
const;
1775 template <
class OutputVector,
typename number>
1777 distribute_local_to_global_by_interpolation(
1779 OutputVector &values,
1796 template <
typename number,
typename OutputVector>
1799 OutputVector &global_destination)
const;
1815 template <
typename ForwardIterator,
typename OutputVector>
1817 distribute_local_to_global(ForwardIterator local_source_begin,
1818 ForwardIterator local_source_end,
1819 OutputVector &global_destination)
const;
1834 template <
typename ForwardIterator,
typename OutputVector>
1836 distribute_local_to_global(
1838 ForwardIterator local_source_begin,
1839 ForwardIterator local_source_end,
1840 OutputVector &global_destination)
const;
1848 template <
typename number,
typename OutputMatrix>
1851 OutputMatrix &global_destination)
const;
1857 template <
typename number,
typename OutputMatrix,
typename OutputVector>
1861 OutputMatrix &global_matrix,
1862 OutputVector &global_vector)
const;
1889 get_active_or_mg_dof_indices(
1890 std::vector<types::global_dof_index> &dof_indices)
const;
1924 get_dof_indices(std::vector<types::global_dof_index> &dof_indices)
const;
1931 get_mg_dof_indices(std::vector<types::global_dof_index> &dof_indices)
const;
1986 active_fe_index()
const;
2025 set_dof_indices(
const std::vector<types::global_dof_index> &dof_indices);
2031 set_mg_dof_indices(
const std::vector<types::global_dof_index> &dof_indices);
2059 get_future_fe()
const;
2080 future_fe_index()
const;
2100 future_fe_index_set()
const;
2111 clear_future_fe_index()
const;
2117 friend struct ::internal::DoFCellAccessorImplementation::Implementation;
2121template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
2125 return level_dof_access;
2130template <
int structdim,
int dim,
int spacedim>
2131template <
typename OtherAccessor>
2133 const OtherAccessor &)
2136 ExcMessage(
"You are attempting an illegal conversion between "
2137 "iterator/accessor types. The constructor you call "
2138 "only exists to make certain template constructs "
2139 "easier to write as dimension independent code but "
2140 "the conversion is not valid in the current context."));
2148template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
2156template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
2162 : ::
internal::DoFAccessorImplementation::
2163 Inheritance<structdim, dim, spacedim>::
BaseClass(tria,
level, index)
2169 "You can't create a DoF accessor in which the DoFHandler object "
2170 "uses a different triangulation than the one you pass as argument."));
2175template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
2176template <
int structdim2,
int dim2,
int spacedim2>
2180 Assert(
false, ExcInvalidObject());
2185template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
2186template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
2190 , dof_handler(nullptr)
2194 "You are trying to assign iterators that are incompatible. "
2195 "The reason for incompatibility is that they refer to objects of "
2196 "different dimensionality (e.g., assigning a line iterator "
2197 "to a quad iterator)."));
2202template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
2203template <
bool level_dof_access2>
2207 , dof_handler(const_cast<
DoFHandler<dim, spacedim> *>(other.dof_handler))
2212template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
2217 Assert(dh !=
nullptr, ExcInvalidObject());
2218 this->dof_handler = dh;
2223template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
2227 Assert(this->dof_handler !=
nullptr, ExcInvalidObject());
2228 return *this->dof_handler;
2233template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
2238 Assert(this->dof_handler !=
nullptr, ExcInvalidObject());
2239 BaseClass::copy_from(da);
2244template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
2245template <
bool level_dof_access2>
2250 BaseClass::copy_from(a);
2256template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
2257template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
2262 Assert(structdim == structdim2, ExcCantCompareIterators());
2264 return (BaseClass::operator==(a));
2269template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
2270template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
2275 Assert(structdim == structdim2, ExcCantCompareIterators());
2277 return (BaseClass::operator!=(a));
2282template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
2285 const unsigned int i)
const
2287 Assert(
static_cast<unsigned int>(this->present_level) <
2288 this->dof_handler->object_dof_indices.size(),
2295 *t, this->dof_handler);
2302 namespace DoFAccessorImplementation
2308 template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
2322 "It is not possible to specify a FE index if no hp support is used!"));
2342 "You need to specify a FE index if hp support is used!"));
2363 boost::container::small_vector<::types::global_dof_index, 27>;
2378 typename GlobalIndexType,
2379 typename DoFPProcessor>
2382 const unsigned int obj_level,
2383 const unsigned int obj_index,
2385 const unsigned int local_index,
2386 const std::integral_constant<int, structdim> &,
2387 GlobalIndexType &global_index,
2388 const DoFPProcessor &process)
2400 [obj_level][structdim]
2409 if (structdim == dim)
2413 [obj_level][structdim]
2441 "You are requesting an active FE index that is not assigned "
2442 "to any of the cells connected to this entity."));
2457 [obj_level][structdim]
2466 [obj_level][structdim]
2468 [obj_level][structdim]
2486 template <
int dim,
int spacedim,
int structdim>
2487 static std::pair<unsigned int, unsigned int>
2489 const unsigned int obj_level,
2490 const unsigned int obj_index,
2492 const std::integral_constant<int, structdim> &)
2498 if (structdim == dim)
2500 const unsigned int ptr_0 =
2502 const unsigned int length =
2503 dof_handler.
get_fe(fe_index).template n_dofs_per_object<dim>(0);
2505 return {ptr_0, length};
2514 const unsigned int ptr_0 =
2516 const unsigned int length =
2517 dof_handler.
object_dof_ptr[obj_level][structdim][obj_index + 1] -
2520 return {ptr_0, length};
2532 const auto fe_index_local_ptr =
2539 Assert(fe_index_local_ptr !=
2543 "You tried to call a function accessing DoF indices, but "
2544 "they appear not be available (yet) or inconsistent. "
2545 "Did you call distribute_dofs() first? Alternatively, if "
2546 "you are using different elements on different cells (i.e., "
2547 "you are using the hp capabilities of deal.II), did you "
2548 "change the active_fe_index of a cell since you last "
2549 "called distribute_dofs()?"));
2554 fe_index_local_ptr);
2560 const unsigned int ptr_0 =
2565 const unsigned int ptr_1 =
2569 fe_index_local + 1];
2571 return {ptr_0, ptr_1 - ptr_0};
2574 template <
int dim,
int spacedim,
int structdim,
bool level_dof_access>
2575 static std::pair<unsigned int, unsigned int>
2577 const ::DoFAccessor<structdim, dim, spacedim, level_dof_access>
2585 std::integral_constant<int, structdim>());
2588 template <
int dim,
int spacedim,
int structdim>
2589 static std::pair<unsigned int, unsigned int>
2611 typename DoFProcessor,
2612 typename DoFMapping>
2615 const unsigned int obj_level,
2616 const unsigned int obj_index,
2618 const DoFMapping &mapping,
2619 const std::integral_constant<int, structdim> &dd,
2621 const DoFProcessor &process)
2628 if (range.second == 0)
2631 std::vector<types::global_dof_index> &object_dof_indices =
2636 object_dof_indices.data() + range.first;
2639 for (
unsigned int i = 0; i < range.second; ++i)
2642 stored_indices[(structdim == 0 || structdim == dim) ? i :
2645 if (dof_indices_ptr !=
nullptr)
2662 template <
int dim,
int spacedim,
int structdim>
2665 const unsigned int obj_level,
2666 const unsigned int obj_index,
2668 const unsigned int local_index,
2669 const std::integral_constant<int, structdim> &dd,
2679 [](
auto &ptr,
const auto &
value) { ptr =
value; });
2693 template <
int dim,
int spacedim,
int structdim>
2696 const unsigned int obj_level,
2697 const unsigned int obj_index,
2699 const unsigned int local_index,
2700 const std::integral_constant<int, structdim> &dd)
2710 [](
const auto &ptr,
auto &
value) {
value = ptr; });
2711 return global_index;
2715 template <
int dim,
int spacedim>
2719 const unsigned int vertex_index,
2720 const unsigned int i)
2724 "DoFHandler in hp-mode does not implement multilevel DoFs."));
2741 template <
int dim,
int spacedim,
int structdim>
2744 const unsigned int obj_level,
2745 const unsigned int obj_index,
2746 const std::integral_constant<int, structdim> &)
2757 if (structdim == dim)
2780 template <
int dim,
int spacedim,
int structdim>
2783 const unsigned int obj_level,
2784 const unsigned int obj_index,
2785 const unsigned int local_index,
2786 const std::integral_constant<int, structdim> &)
2791 Assert(((structdim == dim) &&
2801 if (structdim == dim)
2832 template <
int dim,
int spacedim,
int structdim>
2833 static std::set<types::fe_index>
2835 const unsigned int obj_level,
2836 const unsigned int obj_index,
2837 const std::integral_constant<int, structdim> &t)
2846 if (structdim == dim)
2850 std::set<types::fe_index> active_fe_indices;
2851 for (
unsigned int i = 0;
2854 active_fe_indices.insert(
2856 return active_fe_indices;
2861 template <
int dim,
int spacedim,
int structdim>
2864 const unsigned int obj_level,
2865 const unsigned int obj_index,
2867 const std::integral_constant<int, structdim> &)
2876 if (structdim == dim)
2893 template <
typename InputVector,
typename ForwardIterator>
2898 ForwardIterator local_values_begin)
2900 values.extract_subvector_to(cache, cache_end, local_values_begin);
2905#ifdef DEAL_II_WITH_TRILINOS
2906 static std::vector<unsigned int>
2911 std::vector<unsigned int> idx(v_end - v_begin);
2912 std::iota(idx.begin(), idx.end(), 0u);
2915 std::sort(idx.begin(),
2917 [&v_begin](
unsigned int i1,
unsigned int i2) {
2918 return *(v_begin + i1) < *(v_begin + i2);
2926# ifdef DEAL_II_TRILINOS_WITH_TPETRA
2927 template <
typename ForwardIterator,
typename Number,
typename MemorySpace>
2934 ForwardIterator local_values_begin)
2936 std::vector<unsigned int> sorted_indices_pos =
2938 const unsigned int cache_size = cache_end - cache_begin;
2939 std::vector<types::global_dof_index> cache_indices(cache_size);
2940 for (
unsigned int i = 0; i < cache_size; ++i)
2941 cache_indices[i] = *(cache_begin + sorted_indices_pos[i]);
2943 IndexSet index_set(cache_indices.back() + 1);
2944 index_set.
add_indices(cache_indices.begin(), cache_indices.end());
2950 for (
unsigned int i = 0; i < cache_size; ++i, ++local_values_begin)
2951 *local_values_begin = read_write_vector[sorted_indices_pos[i]];
2957 template <
typename ForwardIterator>
2962 ForwardIterator local_values_begin)
2964 std::vector<unsigned int> sorted_indices_pos =
2966 const unsigned int cache_size = cache_end - cache_begin;
2967 std::vector<types::global_dof_index> cache_indices(cache_size);
2968 for (
unsigned int i = 0; i < cache_size; ++i)
2969 cache_indices[i] = *(cache_begin + sorted_indices_pos[i]);
2971 IndexSet index_set(cache_indices.back() + 1);
2972 index_set.
add_indices(cache_indices.begin(), cache_indices.end());
2978 for (
unsigned int i = 0; i < cache_size; ++i, ++local_values_begin)
2979 *local_values_begin = read_write_vector[sorted_indices_pos[i]];
2987 template <
int dim,
int spacedim,
bool level_dof_access,
int structdim>
2990 const ::DoFAccessor<structdim, dim, spacedim, level_dof_access>
2993 const bool count_level_dofs)
2998 if (count_level_dofs)
3000 const auto &fe = accessor.get_fe(fe_index_);
3003 dofs_per_vertex = fe.n_dofs_per_vertex(),
3004 dofs_per_line = fe.n_dofs_per_line(),
3005 dofs_per_quad = fe.n_dofs_per_quad(0 ),
3006 dofs_per_hex = fe.n_dofs_per_hex();
3008 unsigned int index = 0;
3011 index += dofs_per_vertex * accessor.n_vertices();
3014 if (structdim == 2 || structdim == 3)
3015 index += dofs_per_line * accessor.n_lines();
3019 index += dofs_per_quad * accessor.n_faces();
3022 const unsigned int interior_dofs =
3023 structdim == 1 ? dofs_per_line :
3024 (structdim == 2 ? dofs_per_quad : dofs_per_hex);
3026 index += interior_dofs;
3032 const auto fe_index =
3034 accessor, fe_index_);
3036 unsigned int index = 0;
3039 for (
const auto vertex : accessor.vertex_indices())
3042 accessor.vertex_index(vertex),
3044 std::integral_constant<int, 0>())
3048 if (structdim == 2 || structdim == 3)
3049 for (
const auto line : accessor.line_indices())
3055 for (
const auto face : accessor.face_indices())
3072 template <
typename ArrayType>
3076 return array.size();
3085 template <
typename ArrayType>
3107 bool level_dof_access,
3109 typename DoFIndicesType,
3110 typename DoFOperation,
3111 typename DoFProcessor>
3114 const ::DoFAccessor<structdim, dim, spacedim, level_dof_access>
3116 const DoFIndicesType &const_dof_indices,
3118 const DoFOperation &dof_operation,
3119 const DoFProcessor &dof_processor,
3120 const bool count_level_dofs)
3124 accessor, fe_index_);
3129 (void)count_level_dofs;
3131 const auto &fe = accessor.get_fe(fe_index);
3145 if (fe.n_dofs_per_vertex() > 0)
3146 for (
const auto vertex : accessor.vertex_indices())
3147 dof_operation.process_vertex_dofs(*accessor.dof_handler,
3148 accessor.vertex_index(vertex),
3160 if ((structdim == 2 || structdim == 3) && fe.n_dofs_per_line() > 0)
3162 const auto line_indices = internal::TriaAccessorImplementation::
3163 Implementation::get_line_indices_of_cell(accessor);
3164 const auto line_orientations =
3165 internal::TriaAccessorImplementation::Implementation::
3166 get_line_orientations_of_cell(accessor);
3168 for (
const auto line : accessor.line_indices())
3170 const auto line_orientation = line_orientations[line];
3172 dof_operation.process_dofs(
3173 accessor.get_dof_handler(),
3177 [](
const auto d) {
return d; },
3178 std::integral_constant<int, 1>(),
3182 dof_operation.process_dofs(
3183 accessor.get_dof_handler(),
3187 [&fe, line_orientation](
const auto d) {
3188 return fe.adjust_line_dof_index_for_line_orientation(
3189 d, line_orientation);
3191 std::integral_constant<int, 1>(),
3205 if (structdim == 3 && fe.max_dofs_per_quad() > 0)
3206 for (
const auto face_no : accessor.face_indices())
3208 const auto combined_orientation =
3209 accessor.combined_face_orientation(face_no);
3210 const unsigned int quad_index = accessor.quad_index(face_no);
3211 if (combined_orientation ==
3213 dof_operation.process_dofs(
3214 accessor.get_dof_handler(),
3218 [](
const auto d) {
return d; },
3219 std::integral_constant<int, 2>(),
3223 dof_operation.process_dofs(
3224 accessor.get_dof_handler(),
3229 return fe.adjust_quad_dof_index_for_face_orientation(
3230 d, face_no, combined_orientation);
3232 std::integral_constant<int, 2>(),
3241 if (((dim == 3 && structdim == 2) ?
3242 fe.max_dofs_per_quad() :
3243 fe.template n_dofs_per_object<structdim>()) > 0)
3244 dof_operation.process_dofs(
3245 accessor.get_dof_handler(),
3249 [&](
const auto d) {
return d; },
3250 std::integral_constant<int, structdim>(),
3254 if (dof_indices_ptr !=
nullptr)
3267 for (; dof_indices_ptr < end_dof_indices; ++dof_indices_ptr)
3268 dof_processor(invalid_index, dof_indices_ptr);
3277 template <
int dim,
int spacedim>
3283 template <
typename DoFProcessor>
3286 const unsigned int vertex_index,
3289 const DoFProcessor &dof_processor)
const
3300 std::integral_constant<int, 0>(),
3308 template <
int structdim,
typename DoFMapping,
typename DoFProcessor>
3311 const unsigned int obj_level,
3312 const unsigned int obj_index,
3314 const DoFMapping &mapping,
3315 const std::integral_constant<int, structdim>,
3317 const DoFProcessor &dof_processor)
const
3325 std::integral_constant<
int,
std::min(structdim, dim)>(),
3337 template <
int dim,
int spacedim>
3350 template <
typename DoFProcessor>
3353 const unsigned int vertex_index,
3356 const DoFProcessor &dof_processor)
const
3358 const unsigned int n_indices =
3359 dof_handler.
get_fe(0).template n_dofs_per_object<0>();
3364 for (
unsigned int d = 0; d < n_indices; ++d, ++dof_indices_ptr)
3365 dof_processor(stored_indices[d], dof_indices_ptr);
3371 template <
int structdim,
typename DoFMapping,
typename DoFProcessor>
3375 const unsigned int obj_index,
3377 const DoFMapping &mapping,
3378 const std::integral_constant<int, structdim>,
3380 const DoFProcessor &dof_processor)
const
3382 const unsigned int n_indices =
3383 dof_handler.
get_fe(0).template n_dofs_per_object<structdim>();
3391 std::integral_constant<
int,
std::min(structdim, dim)>());
3392 for (
unsigned int d = 0; d < n_indices; ++d, ++dof_indices_ptr)
3393 dof_processor(stored_indices[structdim < dim ? mapping(d) : d],
3403 template <
int dim,
int spacedim,
bool level_dof_access,
int structdim>
3406 const ::DoFAccessor<structdim, dim, spacedim, level_dof_access>
3408 std::vector<types::global_dof_index> &dof_indices,
3416 [](
auto stored_index,
auto dof_ptr) { *dof_ptr = stored_index; },
3422 template <
int dim,
int spacedim,
bool level_dof_access,
int structdim>
3425 const ::DoFAccessor<structdim, dim, spacedim, level_dof_access>
3427 const std::vector<types::global_dof_index> &dof_indices,
3438 "This function is intended to be used for DoFCellAccessor, i.e., "
3439 "dimension == structdim."));
3446 [](
auto &stored_index,
auto dof_ptr) { stored_index = *dof_ptr; },
3452 template <
int dim,
int spacedim,
bool level_dof_access,
int structdim>
3455 const ::DoFAccessor<structdim, dim, spacedim, level_dof_access>
3458 std::vector<types::global_dof_index> &dof_indices,
3462 ExcMessage(
"MG DoF indices cannot be queried in hp case"));
3468 [](
auto stored_index,
auto dof_ptr) { *dof_ptr = stored_index; },
3474 template <
int dim,
int spacedim,
bool level_dof_access,
int structdim>
3477 const ::DoFAccessor<structdim, dim, spacedim, level_dof_access>
3480 const std::vector<types::global_dof_index> &dof_indices,
3484 ExcMessage(
"MG DoF indices cannot be queried in hp case"));
3492 ExcMessage(
"This function is intended to be used for "
3493 "DoFCellAccessor, i.e., dimension == structdim."));
3500 [](
auto &stored_index,
auto dof_ptr) { stored_index = *dof_ptr; },
3506 template <
int dim,
int spacedim>
3514 const unsigned int obj_index,
3516 const unsigned int local_index,
3517 const std::integral_constant<int, dim>)
3522 return mg_level->dof_object.access_dof_index(
3531 template <
int dim,
int spacedim, std::enable_if_t<(dim > 1),
int> = 0>
3539 const unsigned int obj_index,
3541 const unsigned int local_index,
3542 const std::integral_constant<int, 1>)
3544 return mg_faces->lines.access_dof_index(
3553 template <
int spacedim>
3561 const unsigned int obj_index,
3563 const unsigned int local_index,
3564 const std::integral_constant<int, 2>)
3568 return mg_faces->quads.access_dof_index(
3578 template <
int dim,
int spacedim,
bool level_dof_access>
3581 const ::DoFCellAccessor<dim, spacedim, level_dof_access> &accessor,
3583 const unsigned int fe_index);
3589template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
3592 const unsigned int i,
3595 const auto fe_index =
3600 return ::internal::DoFAccessorImplementation::Implementation::
3601 get_dof_index(*this->dof_handler,
3606 std::integral_constant<int, structdim>());
3610template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
3614 const unsigned int i)
const
3618 this->dof_handler->mg_levels[
level],
3619 this->dof_handler->mg_faces,
3623 std::integral_constant<int, structdim>());
3627template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
3630 const unsigned int i,
3634 const auto fe_index =
3645 std::integral_constant<int, structdim>(),
3651template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
3655 const unsigned int i,
3660 this->dof_handler->mg_levels[
level],
3661 this->dof_handler->mg_faces,
3665 std::integral_constant<int, structdim>()) = index;
3670template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
3676 return ::internal::DoFAccessorImplementation::Implementation::
3677 n_active_fe_indices(*this->dof_handler,
3680 std::integral_constant<int, structdim>());
3685template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
3688 const unsigned int n)
const
3691 return ::internal::DoFAccessorImplementation::Implementation::
3692 nth_active_fe_index(*this->dof_handler,
3696 std::integral_constant<int, structdim>());
3701template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
3702inline std::set<types::fe_index>
3706 std::set<types::fe_index> active_fe_indices;
3707 for (
unsigned int i = 0; i < n_active_fe_indices(); ++i)
3708 active_fe_indices.insert(nth_active_fe_index(i));
3709 return active_fe_indices;
3714template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
3720 return ::internal::DoFAccessorImplementation::Implementation::
3721 fe_index_is_active(*this->dof_handler,
3725 std::integral_constant<int, structdim>());
3730template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
3733 const unsigned int vertex,
3734 const unsigned int i,
3738 (((this->dof_handler->hp_capability_enabled ==
false) &&
3749 this->nth_active_fe_index(0) :
3752 return ::internal::DoFAccessorImplementation::Implementation::
3753 get_dof_index(*this->dof_handler,
3755 this->vertex_index(vertex),
3758 std::integral_constant<int, 0>());
3762template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
3766 const unsigned int vertex,
3767 const unsigned int i,
3770 const auto fe_index =
3774 Assert(this->dof_handler !=
nullptr, ExcInvalidObject());
3775 Assert(this->dof_handler->mg_vertex_dofs.size() > 0,
3776 ExcMessage(
"Multigrid DoF indices can only be accessed after "
3777 "DoFHandler::distribute_mg_dofs() has been called!"));
3779 AssertIndexRange(i, this->dof_handler->get_fe(fe_index).n_dofs_per_vertex());
3781 Assert(dof_handler->hp_capability_enabled ==
false,
3783 "DoFHandler in hp-mode does not implement multilevel DoFs."));
3785 return this->dof_handler->mg_vertex_dofs[this->vertex_index(vertex)]
3786 .access_index(
level, i, this->dof_handler->get_fe().n_dofs_per_vertex());
3791template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
3795 const unsigned int vertex,
3796 const unsigned int i,
3800 const auto fe_index =
3804 Assert(this->dof_handler !=
nullptr, ExcInvalidObject());
3806 AssertIndexRange(i, this->dof_handler->get_fe(fe_index).n_dofs_per_vertex());
3808 Assert(dof_handler->hp_capability_enabled ==
false,
3810 "DoFHandler in hp-mode does not implement multilevel DoFs."));
3812 this->dof_handler->mg_vertex_dofs[this->vertex_index(vertex)].access_index(
3813 level, i, this->dof_handler->get_fe().n_dofs_per_vertex()) = index;
3818template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
3823 Assert(fe_index_is_active(fe_index) ==
true,
3824 ExcMessage(
"This function can only be called for active FE indices"));
3826 return this->dof_handler->get_fe(fe_index);
3831template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
3832inline typename ::internal::DoFHandlerImplementation::
3833 Iterators<dim, spacedim, level_dof_access>::line_iterator
3835 const unsigned int i)
const
3843 ExcMessage(
"You can only ask for line zero if the "
3844 "current object is a line itself."));
3845 return typename ::internal::DoFHandlerImplementation::
3846 Iterators<dim, spacedim, level_dof_access>::cell_iterator(
3847 &this->get_triangulation(),
3850 &this->get_dof_handler());
3858 return typename ::internal::DoFHandlerImplementation::
3859 Iterators<dim, spacedim, level_dof_access>::line_iterator(
3862 this->line_index(i),
3867template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
3868inline typename ::internal::DoFHandlerImplementation::
3869 Iterators<dim, spacedim, level_dof_access>::quad_iterator
3871 const unsigned int i)
const
3881 ExcMessage(
"You can only ask for quad zero if the "
3882 "current object is a quad itself."));
3883 return typename ::internal::DoFHandlerImplementation::
3884 Iterators<dim, spacedim>::cell_iterator(&this->get_triangulation(),
3887 &this->get_dof_handler());
3895 return typename ::internal::DoFHandlerImplementation::
3896 Iterators<dim, spacedim, level_dof_access>::quad_iterator(
3899 this->quad_index(i),
3907template <
int spacedim,
bool level_dof_access>
3910 Assert(
false, ExcInvalidObject());
3915template <
int spacedim,
bool level_dof_access>
3919 const unsigned int vertex_index,
3921 :
BaseClass(tria, vertex_kind, vertex_index)
3922 , dof_handler(const_cast<
DoFHandler<1, spacedim> *>(dof_handler))
3927template <
int spacedim,
bool level_dof_access>
3949 Assert((tria ==
nullptr) && (
level == -2) && (index == -2) &&
3952 "This constructor can not be called for face iterators in 1d, "
3953 "except to default-construct iterator objects."));
3958template <
int spacedim,
bool level_dof_access>
3959template <
int structdim2,
int dim2,
int spacedim2>
3963 Assert(
false, ExcInvalidObject());
3968template <
int spacedim,
bool level_dof_access>
3969template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
3973 Assert(
false, ExcInvalidObject());
3978template <
int spacedim,
bool level_dof_access>
3983 Assert(dh !=
nullptr, ExcInvalidObject());
3984 this->dof_handler = dh;
3989template <
int spacedim,
bool level_dof_access>
3992 const unsigned int ,
4001template <
int spacedim,
bool level_dof_access>
4005 return *this->dof_handler;
4010template <
int spacedim,
bool level_dof_access>
4013 std::vector<types::global_dof_index> &dof_indices,
4016 const auto fe_index =
4020 for (
unsigned int i = 0; i < dof_indices.size(); ++i)
4021 dof_indices[i] = ::internal::DoFAccessorImplementation::
4022 Implementation::get_dof_index(*dof_handler,
4024 this->global_vertex_index,
4027 std::integral_constant<int, 0>());
4032template <
int spacedim,
bool level_dof_access>
4036 std::vector<types::global_dof_index> &dof_indices,
4039 const auto fe_index =
4045 for (
unsigned int i = 0; i < dof_indices.size(); ++i)
4048 mg_vertex_dof_index(*dof_handler,
level, this->global_vertex_index, i);
4053template <
int spacedim,
bool level_dof_access>
4056 const unsigned int vertex,
4057 const unsigned int i,
4060 const auto fe_index =
4066 return ::internal::DoFAccessorImplementation::Implementation::
4067 get_dof_index(*dof_handler,
4069 this->global_vertex_index,
4072 std::integral_constant<int, 0>());
4077template <
int spacedim,
bool level_dof_access>
4080 const unsigned int i,
4083 const auto fe_index =
4087 return ::internal::DoFAccessorImplementation::Implementation::
4088 get_dof_index(*this->dof_handler,
4090 this->vertex_index(0),
4093 std::integral_constant<int, 0>());
4098template <
int spacedim,
bool level_dof_access>
4107template <
int spacedim,
bool level_dof_access>
4110 const unsigned int )
const
4117template <
int spacedim,
bool level_dof_access>
4128template <
int spacedim,
bool level_dof_access>
4133 Assert(this->dof_handler !=
nullptr, ExcInvalidObject());
4134 return dof_handler->get_fe(fe_index);
4139template <
int spacedim,
bool level_dof_access>
4144 Assert(this->dof_handler !=
nullptr, ExcInvalidObject());
4145 BaseClass::copy_from(da);
4150template <
int spacedim,
bool level_dof_access>
4151template <
bool level_dof_access2>
4156 BaseClass::copy_from(a);
4162template <
int spacedim,
bool level_dof_access>
4165 const unsigned int )
const
4172template <
int spacedim,
bool level_dof_access>
4173inline typename ::internal::DoFHandlerImplementation::
4174 Iterators<1, spacedim, level_dof_access>::line_iterator
4176 const unsigned int )
const
4179 return typename ::internal::DoFHandlerImplementation::
4180 Iterators<1, spacedim, level_dof_access>::line_iterator();
4185template <
int spacedim,
bool level_dof_access>
4186inline typename ::internal::DoFHandlerImplementation::
4187 Iterators<1, spacedim, level_dof_access>::quad_iterator
4189 const unsigned int )
const
4192 return typename ::internal::DoFHandlerImplementation::
4193 Iterators<1, spacedim, level_dof_access>::quad_iterator();
4198template <
int spacedim,
bool level_dof_access>
4199template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
4204 Assert(structdim2 == 0, ExcCantCompareIterators());
4206 return (BaseClass::operator==(a));
4211template <
int spacedim,
bool level_dof_access>
4212template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
4217 Assert(structdim2 == 0, ExcCantCompareIterators());
4219 return (BaseClass::operator!=(a));
4229 namespace DoFCellAccessorImplementation
4241 template <
int dim,
int spacedim,
bool level_dof_access>
4250 (
typename std::decay_t<
decltype(accessor)>::ExcInvalidObject()));
4251 Assert(
static_cast<unsigned int>(accessor.level()) <
4265 template <
int dim,
int spacedim,
bool level_dof_access>
4278 (
typename std::decay_t<
decltype(accessor)>::ExcInvalidObject()));
4279 Assert(
static_cast<unsigned int>(accessor.level()) <
4283 ExcMessage(
"Invalid finite element index."));
4287 [accessor.present_index] = i;
4296 template <
int dim,
int spacedim,
bool level_dof_access>
4305 (
typename std::decay_t<
decltype(accessor)>::ExcInvalidObject()));
4306 Assert(
static_cast<unsigned int>(accessor.level()) <
4313 [accessor.present_index];
4317 [accessor.present_index];
4325 template <
int dim,
int spacedim,
bool level_dof_access>
4338 (
typename std::decay_t<
decltype(accessor)>::ExcInvalidObject()));
4339 Assert(
static_cast<unsigned int>(accessor.level()) <
4343 ExcMessage(
"Invalid finite element index."));
4347 [accessor.present_index] = i;
4356 template <
int dim,
int spacedim,
bool level_dof_access>
4365 (
typename std::decay_t<
decltype(accessor)>::ExcInvalidObject()));
4366 Assert(
static_cast<unsigned int>(accessor.level()) <
4372 [accessor.present_index] !=
4382 template <
int dim,
int spacedim,
bool level_dof_access>
4391 (
typename std::decay_t<
decltype(accessor)>::ExcInvalidObject()));
4392 Assert(
static_cast<unsigned int>(accessor.level()) <
4398 [accessor.present_index] =
4229 namespace DoFCellAccessorImplementation {
…}
4407template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4413 :
DoFAccessor<dimension_, dimension_, space_dimension_, level_dof_access>(
4422template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4423template <
int structdim2,
int dim2,
int spacedim2>
4427 Assert(
false,
typename BaseClass::ExcInvalidObject());
4432template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4433template <
int structdim2,
int dim2,
int spacedim2,
bool level_dof_access2>
4442template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4446 const unsigned int i)
const
4450 this->neighbor_level(i),
4451 this->neighbor_index(i),
4464template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4468 const unsigned int i)
const
4471 q(this->tria, this->
level() + 1, this->child_index(i), this->dof_handler);
4483template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4484inline boost::container::small_vector<
4490 boost::container::small_vector<
4494 child_iterators(this->n_children());
4496 for (
unsigned int i = 0; i < this->n_children(); ++i)
4497 child_iterators[i] = this->child(i);
4499 return child_iterators;
4504template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4510 q(this->tria, this->
level() - 1, this->parent_index(), this->dof_handler);
4519 namespace DoFCellAccessorImplementation
4521 template <
int dim,
int spacedim,
bool level_dof_access>
4525 const ::DoFCellAccessor<dim, spacedim, level_dof_access> &cell,
4526 const unsigned int i,
4527 const std::integral_constant<int, 1>)
4530 &cell.get_triangulation(),
4531 ((i == 0) && cell.at_boundary(0) ?
4533 ((i == 1) && cell.at_boundary(1) ?
4536 cell.vertex_index(i),
4537 &cell.get_dof_handler());
4538 return ::TriaIterator<
4543 template <
int dim,
int spacedim,
bool level_dof_access>
4547 const ::DoFCellAccessor<dim, spacedim, level_dof_access> &cell,
4548 const unsigned int i,
4549 const std::integral_constant<int, 2>)
4551 return cell.
line(i);
4555 template <
int dim,
int spacedim,
bool level_dof_access>
4559 const ::DoFCellAccessor<dim, spacedim, level_dof_access> &cell,
4560 const unsigned int i,
4561 const std::integral_constant<int, 3>)
4563 return cell.
quad(i);
4570template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4573 level_dof_access>::face_iterator
4575 const unsigned int i)
const
4579 return ::internal::DoFCellAccessorImplementation::get_face(
4580 *
this, i, std::integral_constant<int, dimension_>());
4585template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4586inline boost::container::small_vector<
4593 boost::container::small_vector<
4597 face_iterators(this->n_faces());
4599 for (
const unsigned int i : this->face_indices())
4602 *
this, i, std::integral_constant<int, dimension_>());
4604 return face_iterators;
4609template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4613 std::vector<types::global_dof_index> &dof_indices)
const
4615 if (level_dof_access)
4616 get_mg_dof_indices(dof_indices);
4618 get_dof_indices(dof_indices);
4623template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4624template <
class InputVector,
typename number>
4627 const InputVector &values,
4630 get_dof_values(values, local_values.
begin(), local_values.
end());
4635template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4636template <
typename Number,
typename ForwardIterator>
4640 ForwardIterator local_values_begin,
4641 ForwardIterator local_values_end)
const
4643 (void)local_values_end;
4644 Assert(this->is_artificial() ==
false,
4645 ExcMessage(
"Can't ask for DoF indices on artificial cells."));
4647 Assert(this->dof_handler !=
nullptr,
typename BaseClass::ExcInvalidObject());
4649 Assert(
static_cast<unsigned int>(local_values_end - local_values_begin) ==
4650 this->get_fe().n_dofs_per_cell(),
4652 Assert(values.size() == this->get_dof_handler().n_dofs(),
4656 dof_indices(this->get_fe().n_dofs_per_cell());
4658 *
this, dof_indices, this->active_fe_index());
4660 boost::container::small_vector<Number, 27> values_temp(local_values_end -
4661 local_values_begin);
4666 using view_type = std::remove_reference_t<
decltype(*local_values_begin)>;
4668 local_values_end - local_values_begin);
4669 std::copy(values_temp.begin(), values_temp.end(), values_view2.
begin());
4674template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4675template <
class InputVector,
typename ForwardIterator>
4679 const InputVector &values,
4680 ForwardIterator local_values_begin,
4681 ForwardIterator local_values_end)
const
4683 Assert(this->is_artificial() ==
false,
4684 ExcMessage(
"Can't ask for DoF indices on artificial cells."));
4687 Assert(
static_cast<unsigned int>(local_values_end - local_values_begin) ==
4688 this->get_fe().n_dofs_per_cell(),
4690 Assert(values.size() == this->get_dof_handler().n_dofs(),
4695 dof_indices(this->get_fe().n_dofs_per_cell());
4697 *
this, dof_indices, this->active_fe_index());
4707template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4708template <
class OutputVector,
typename number>
4712 OutputVector &values)
const
4714 Assert(this->is_artificial() ==
false,
4715 ExcMessage(
"Can't ask for DoF indices on artificial cells."));
4718 Assert(
static_cast<unsigned int>(local_values.
size()) ==
4719 this->get_fe().n_dofs_per_cell(),
4721 Assert(values.size() == this->get_dof_handler().n_dofs(),
4725 Assert(this->dof_handler !=
nullptr,
typename BaseClass::ExcInvalidObject());
4727 dof_indices(this->get_fe().n_dofs_per_cell());
4729 *
this, dof_indices, this->active_fe_index());
4731 for (
unsigned int i = 0; i < this->get_fe().n_dofs_per_cell(); ++i)
4739template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4743 Assert(this->dof_handler !=
nullptr,
typename BaseClass::ExcInvalidObject());
4744 Assert((this->dof_handler->hp_capability_enabled ==
false) ||
4747 "For DoFHandler objects in hp-mode, finite elements are only "
4748 "associated with active cells. Consequently, you can not ask "
4749 "for the active finite element on cells with children."));
4751 const auto &fe = this->dof_handler->get_fe(active_fe_index());
4753 Assert(this->reference_cell() == fe.reference_cell(),
4755 fe.reference_cell()));
4762template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4767 Assert((this->dof_handler->hp_capability_enabled ==
false) ||
4770 "You can not ask for the active FE index on a cell that has "
4771 "children because no degrees of freedom are assigned "
4772 "to this cell and, consequently, no finite element "
4773 "is associated with it."));
4774 Assert((this->dof_handler->hp_capability_enabled ==
false) ||
4775 (this->is_locally_owned() || this->is_ghost()),
4776 ExcMessage(
"You can only query active FE index information on cells "
4777 "that are either locally owned or (after distributing "
4778 "degrees of freedom) are ghost cells."));
4780 return ::internal::DoFCellAccessorImplementation::Implementation::
4781 active_fe_index(*
this);
4786template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4791 Assert((this->dof_handler->hp_capability_enabled ==
false) ||
4793 ExcMessage(
"You can not set the active FE index on a cell that has "
4794 "children because no degrees of freedom will be assigned "
4797 Assert((this->dof_handler->hp_capability_enabled ==
false) ||
4798 this->is_locally_owned(),
4799 ExcMessage(
"You can only set active FE index information on cells "
4800 "that are locally owned. On ghost cells, this information "
4801 "will automatically be propagated from the owning process "
4802 "of that cell, and there is no information at all on "
4803 "artificial cells."));
4811template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4816 Assert(this->dof_handler !=
nullptr,
typename BaseClass::ExcInvalidObject());
4817 Assert((this->dof_handler->hp_capability_enabled ==
false) ||
4820 "For DoFHandler objects in hp-mode, finite elements are only "
4821 "associated with active cells. Consequently, you can not ask "
4822 "for the future finite element on cells with children."));
4824 return this->dof_handler->get_fe(future_fe_index());
4829template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4834 Assert((this->dof_handler->hp_capability_enabled ==
false) ||
4835 (this->has_children() ==
false),
4837 "You can not ask for the future FE index on a cell that has "
4838 "children because no degrees of freedom are assigned "
4839 "to this cell and, consequently, no finite element "
4840 "is associated with it."));
4841 Assert((this->dof_handler->hp_capability_enabled ==
false) ||
4842 (this->is_locally_owned()),
4843 ExcMessage(
"You can only query future FE index information on cells "
4844 "that are locally owned."));
4846 return ::internal::DoFCellAccessorImplementation::Implementation::
4847 future_fe_index(*
this);
4852template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4857 Assert((this->dof_handler->hp_capability_enabled ==
false) ||
4858 (this->has_children() ==
false),
4859 ExcMessage(
"You can not set the future FE index on a cell that has "
4860 "children because no degrees of freedom will be assigned "
4863 Assert((this->dof_handler->hp_capability_enabled ==
false) ||
4864 this->is_locally_owned(),
4865 ExcMessage(
"You can only set future FE index information on cells "
4866 "that are locally owned."));
4874template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4879 Assert((this->dof_handler->hp_capability_enabled ==
false) ||
4880 (this->has_children() ==
false),
4882 "You can not ask for the future FE index on a cell that has "
4883 "children because no degrees of freedom are assigned "
4884 "to this cell and, consequently, no finite element "
4885 "is associated with it."));
4886 Assert((this->dof_handler->hp_capability_enabled ==
false) ||
4887 (this->is_locally_owned()),
4888 ExcMessage(
"You can only query future FE index information on cells "
4889 "that are locally owned."));
4891 return ::internal::DoFCellAccessorImplementation::Implementation::
4892 future_fe_index_set(*
this);
4897template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4902 Assert((this->dof_handler->hp_capability_enabled ==
false) ||
4903 (this->has_children() ==
false),
4905 "You can not ask for the future FE index on a cell that has "
4906 "children because no degrees of freedom are assigned "
4907 "to this cell and, consequently, no finite element "
4908 "is associated with it."));
4909 Assert((this->dof_handler->hp_capability_enabled ==
false) ||
4910 (this->is_locally_owned()),
4911 ExcMessage(
"You can only query future FE index information on cells "
4912 "that are locally owned."));
4920template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4921template <
typename number,
typename OutputVector>
4925 OutputVector &global_destination)
const
4927 this->distribute_local_to_global(local_source.
begin(),
4929 global_destination);
4934template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4935template <
typename ForwardIterator,
typename OutputVector>
4939 ForwardIterator local_source_end,
4940 OutputVector &global_destination)
const
4942 Assert(this->dof_handler !=
nullptr,
4943 (
typename std::decay_t<
decltype(*
this)>::ExcInvalidObject()));
4944 Assert(
static_cast<unsigned int>(local_source_end - local_source_begin) ==
4945 this->get_fe().n_dofs_per_cell(),
4946 (
typename std::decay_t<
decltype(*
this)>::ExcVectorDoesNotMatch()));
4947 Assert(this->dof_handler->n_dofs() == global_destination.size(),
4948 (
typename std::decay_t<
decltype(*
this)>::ExcVectorDoesNotMatch()));
4953 internal::ArrayViewHelper::is_contiguous(local_source_begin,
4956 "This function can not be called with iterator types that do not point to contiguous memory."));
4958 const unsigned int n_dofs = local_source_end - local_source_begin;
4961 dof_indices(n_dofs);
4963 *
this, dof_indices, this->active_fe_index());
4966 global_destination.add(n_dofs, dof_indices.data(), &(*local_source_begin));
4971template <
int dimension_,
int space_dimension_,
bool level_dof_access>
4972template <
typename ForwardIterator,
typename OutputVector>
4977 ForwardIterator local_source_begin,
4978 ForwardIterator local_source_end,
4979 OutputVector &global_destination)
const
4981 Assert(this->dof_handler !=
nullptr,
4982 (
typename std::decay_t<
decltype(*
this)>::ExcInvalidObject()));
4983 Assert(local_source_end - local_source_begin ==
4984 this->get_fe().n_dofs_per_cell(),
4985 (
typename std::decay_t<
decltype(*
this)>::ExcVectorDoesNotMatch()));
4986 Assert(this->dof_handler->n_dofs() == global_destination.size(),
4987 (
typename std::decay_t<
decltype(*
this)>::ExcVectorDoesNotMatch()));
4992 dof_indices(this->get_fe().n_dofs_per_cell());
4994 *
this, dof_indices, this->active_fe_index());
5000 global_destination);
5005template <
int dimension_,
int space_dimension_,
bool level_dof_access>
5006template <
typename number,
typename OutputMatrix>
5010 OutputMatrix &global_destination)
const
5012 Assert(this->dof_handler !=
nullptr,
5013 (
typename std::decay_t<
decltype(*
this)>::ExcInvalidObject()));
5014 Assert(local_source.
m() == this->get_fe().n_dofs_per_cell(),
5015 (
typename std::decay_t<
decltype(*
this)>::ExcMatrixDoesNotMatch()));
5016 Assert(local_source.
n() == this->get_fe().n_dofs_per_cell(),
5017 (
typename std::decay_t<
decltype(*
this)>::ExcMatrixDoesNotMatch()));
5018 Assert(this->dof_handler->n_dofs() == global_destination.m(),
5019 (
typename std::decay_t<
decltype(*
this)>::ExcMatrixDoesNotMatch()));
5020 Assert(this->dof_handler->n_dofs() == global_destination.n(),
5021 (
typename std::decay_t<
decltype(*
this)>::ExcMatrixDoesNotMatch()));
5025 const unsigned int n_dofs = local_source.
m();
5028 dof_indices(n_dofs);
5030 *
this, dof_indices, this->active_fe_index());
5033 for (
unsigned int i = 0; i < n_dofs; ++i)
5034 global_destination.add(dof_indices[i],
5037 &local_source(i, 0));
5042template <
int dimension_,
int space_dimension_,
bool level_dof_access>
5043template <
typename number,
typename OutputMatrix,
typename OutputVector>
5048 OutputMatrix &global_matrix,
5049 OutputVector &global_vector)
const
5051 Assert(this->dof_handler !=
nullptr,
5052 (
typename std::decay_t<
decltype(*
this)>::ExcInvalidObject()));
5053 Assert(local_matrix.
m() == this->get_fe().n_dofs_per_cell(),
5054 (
typename std::decay_t<
decltype(*
this)>::ExcMatrixDoesNotMatch()));
5055 Assert(local_matrix.
n() == this->get_fe().n_dofs_per_cell(),
5056 (
typename std::decay_t<
decltype(*
this)>::ExcVectorDoesNotMatch()));
5057 Assert(this->dof_handler->n_dofs() == global_matrix.m(),
5058 (
typename std::decay_t<
decltype(*
this)>::ExcMatrixDoesNotMatch()));
5059 Assert(this->dof_handler->n_dofs() == global_matrix.n(),
5060 (
typename std::decay_t<
decltype(*
this)>::ExcMatrixDoesNotMatch()));
5061 Assert(local_vector.
size() == this->get_fe().n_dofs_per_cell(),
5062 (
typename std::decay_t<
decltype(*
this)>::ExcVectorDoesNotMatch()));
5063 Assert(this->dof_handler->n_dofs() == global_vector.size(),
5064 (
typename std::decay_t<
decltype(*
this)>::ExcVectorDoesNotMatch()));
5068 const unsigned int n_dofs = this->get_fe().n_dofs_per_cell();
5070 dof_indices(n_dofs);
5072 *
this, dof_indices, this->active_fe_index());
5075 for (
unsigned int i = 0; i < n_dofs; ++i)
5077 global_matrix.add(dof_indices[i],
5080 &local_matrix(i, 0));
5081 global_vector(dof_indices[i]) += local_vector(i);
401 quad(
const unsigned int i)
const; {
…}
392 line(
const unsigned int i)
const; {
…}
ArrayView< std::remove_reference_t< typename std::iterator_traits< Iterator >::reference >, MemorySpaceType > make_array_view(const Iterator begin, const Iterator end)
void get_dof_values(const VectorType &global_vector, ForwardIteratorInd local_indices_begin, ForwardIteratorVec local_vector_begin, ForwardIteratorVec local_vector_end) const
void distribute_local_to_global(const InVector &local_vector, const std::vector< size_type > &local_dof_indices, OutVector &global_vector) const
DoFAccessor< 0, 1, spacedim, level_dof_access > & operator=(const DoFAccessor< 0, 1, spacedim, level_dof_access > &da)=delete
DoFAccessor< 0, 1, spacedim, level_dof_access > & operator=(DoFAccessor< 0, 1, spacedim, level_dof_access > &&) noexcept=default
DoFAccessor(DoFAccessor< 0, 1, spacedim, level_dof_access > &&)=default
DoFAccessor(const DoFAccessor< 0, 1, spacedim, level_dof_access > &)=default
typename::internal::DoFHandlerImplementation::Iterators< dim, spacedim, level_dof_access >::line_iterator line(const unsigned int i) const
friend struct ::internal::DoFHandlerImplementation::Policy::Implementation
static constexpr unsigned int space_dimension
const DoFHandler< dim, spacedim > & get_dof_handler() const
void set_mg_dof_index(const int level, const unsigned int i, const types::global_dof_index index) const
DoFAccessor< structdim, dim, spacedim, level_dof_access > & operator=(const DoFAccessor< structdim, dim, spacedim, level_dof_access > &da)=delete
types::fe_index nth_active_fe_index(const unsigned int n) const
DoFAccessor(DoFAccessor< structdim, dim, spacedim, level_dof_access > &&)
DoFAccessor(const DoFAccessor< structdim, dim, spacedim, level_dof_access2 > &)
TriaIterator< DoFAccessor< structdim, dim, spacedim, level_dof_access > > child(const unsigned int c) const
void set_mg_vertex_dof_index(const int level, const unsigned int vertex, const unsigned int i, const types::global_dof_index index, const types::fe_index fe_index=numbers::invalid_fe_index) const
static constexpr unsigned int dimension
void set_dof_handler(DoFHandler< dim, spacedim > *dh)
DoFAccessor(const DoFAccessor< structdim2, dim2, spacedim2, level_dof_access2 > &)
DoFAccessor< structdim, dim, spacedim, level_dof_access > & operator=(DoFAccessor< structdim, dim, spacedim, level_dof_access > &&)
std::set< types::fe_index > get_active_fe_indices() const
DoFAccessor(const DoFAccessor< structdim, dim, spacedim, level_dof_access > &)=default
void set_mg_dof_indices(const int level, const std::vector< types::global_dof_index > &dof_indices, const types::fe_index fe_index=numbers::invalid_fe_index)
void get_mg_dof_indices(const int level, std::vector< types::global_dof_index > &dof_indices, const types::fe_index fe_index=numbers::invalid_fe_index) const
types::global_dof_index mg_vertex_dof_index(const int level, const unsigned int vertex, const unsigned int i, const types::fe_index fe_index=numbers::invalid_fe_index) const
DoFAccessor(const Triangulation< dim, spacedim > *tria, const int level, const int index, const DoFHandler< dim, spacedim > *dof_handler)
bool operator!=(const DoFAccessor< structdim2, dim2, spacedim2, level_dof_access2 > &) const
void copy_from(const TriaAccessorBase< structdim, dim, spacedim > &da)
types::global_dof_index dof_index(const unsigned int i, const types::fe_index fe_index=numbers::invalid_fe_index) const
typename::internal::DoFHandlerImplementation::Iterators< dim, spacedim, level_dof_access >::quad_iterator quad(const unsigned int i) const
void copy_from(const DoFAccessor< structdim, dim, spacedim, level_dof_access2 > &a)
const FiniteElement< dim, spacedim > & get_fe(const types::fe_index fe_index) const
void get_dof_indices(std::vector< types::global_dof_index > &dof_indices, const types::fe_index fe_index=numbers::invalid_fe_index) const
bool operator==(const DoFAccessor< structdim2, dim2, spacedim2, level_dof_access2 > &) const
types::global_dof_index mg_dof_index(const int level, const unsigned int i) const
void set_dof_index(const unsigned int i, const types::global_dof_index index, const types::fe_index fe_index=numbers::invalid_fe_index) const
unsigned int n_active_fe_indices() const
static bool is_level_cell()
typename ::internal::DoFAccessorImplementation::Inheritance< structdim, dimension, space_dimension >::BaseClass BaseClass
types::global_dof_index vertex_dof_index(const unsigned int vertex, const unsigned int i, const types::fe_index fe_index=numbers::invalid_fe_index) const
DoFHandler< dim, spacedim > * dof_handler
DoFAccessor(const InvalidAccessor< structdim2, dim2, spacedim2 > &)
bool fe_index_is_active(const types::fe_index fe_index) const
boost::container::small_vector< TriaIterator< DoFCellAccessor< dimension_, space_dimension_, level_dof_access > >, GeometryInfo< dimension_ >::max_children_per_cell > child_iterators() const
DoFCellAccessor(DoFCellAccessor< dimension_, space_dimension_, level_dof_access > &&)=default
void get_active_or_mg_dof_indices(std::vector< types::global_dof_index > &dof_indices) const
void get_dof_values(const InputVector &values, Vector< number > &local_values) const
boost::container::small_vector< face_iterator, GeometryInfo< dimension_ >::faces_per_cell > face_iterators() const
types::fe_index future_fe_index() const
DoFCellAccessor< dimension_, space_dimension_, level_dof_access > & operator=(const DoFCellAccessor< dimension_, space_dimension_, level_dof_access > &da)=delete
const FiniteElement< dimension_, space_dimension_ > & get_fe() const
TriaIterator< DoFCellAccessor< dimension_, space_dimension_, level_dof_access > > child(const unsigned int i) const
DoFCellAccessor(const Triangulation< dimension_, space_dimension_ > *tria, const int level, const int index, const AccessorData *local_data)
void set_future_fe_index(const types::fe_index i) const
void distribute_local_to_global(const Vector< number > &local_source, OutputVector &global_destination) const
const FiniteElement< dimension_, space_dimension_ > & get_future_fe() const
TriaIterator< DoFCellAccessor< dimension_, space_dimension_, level_dof_access > > neighbor(const unsigned int i) const
void set_active_fe_index(const types::fe_index i) const
face_iterator face(const unsigned int i) const
void clear_future_fe_index() const
void set_dof_values(const Vector< number > &local_values, OutputVector &values) const
DoFCellAccessor(const DoFCellAccessor< dimension_, space_dimension_, level_dof_access > &)=default
~DoFCellAccessor()=default
bool future_fe_index_set() const
DoFCellAccessor< dimension_, space_dimension_, level_dof_access > & operator=(DoFCellAccessor< dimension_, space_dimension_, level_dof_access > &&)=default
types::fe_index active_fe_index() const
TriaIterator< DoFCellAccessor< dimension_, space_dimension_, level_dof_access > > parent() const
std::vector< std::unique_ptr<::internal::DoFHandlerImplementation::DoFLevel< dim > > > mg_levels
std::vector< std::vector< types::fe_index > > hp_cell_future_fe_indices
const FiniteElement< dim, spacedim > & get_fe(const types::fe_index index=0) const
std::vector< std::vector< types::fe_index > > hp_cell_active_fe_indices
std::vector< MGVertexDoFs > mg_vertex_dofs
std::vector< std::array< std::vector< types::global_dof_index >, dim+1 > > object_dof_indices
const Triangulation< dim, spacedim > & get_triangulation() const
std::vector< std::array< std::vector< offset_type >, dim+1 > > object_dof_ptr
std::unique_ptr<::internal::DoFHandlerImplementation::DoFFaces< dim > > mg_faces
bool hp_capability_enabled
bool has_hp_capabilities() const
std::array< std::vector< offset_type >, dim+1 > hp_object_fe_ptr
std::array< std::vector< types::fe_index >, dim+1 > hp_object_fe_indices
DoFInvalidAccessor(const void *parent=nullptr, const int level=-1, const int index=-1, const AccessorData *local_data=nullptr)
typename InvalidAccessor< structdim, dim, spacedim >::AccessorData AccessorData
unsigned int n_dofs_per_vertex() const
void add_indices(const ForwardIterator &begin, const ForwardIterator &end)
void import_elements(const ::Vector< Number > &vec, VectorOperation::values operation, const std::shared_ptr< const Utilities::MPI::CommunicationPatternBase > &communication_pattern={})
TriaIterator< TriaAccessor< structdim, dim, spacedim > > child(const unsigned int i) const
IteratorState::IteratorStates state() const
virtual size_type size() const override
#define DEAL_II_ALWAYS_INLINE
#define DEAL_II_NAMESPACE_OPEN
constexpr bool running_in_debug_mode()
#define DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_ASSERT_UNREACHABLE()
#define DEAL_II_NOT_IMPLEMENTED()
static ::ExceptionBase & ExcVectorNotEmpty()
#define DeclException0(Exception0)
static ::ExceptionBase & ExcNonMatchingReferenceCellTypes(ReferenceCell arg1, ReferenceCell arg2)
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
static ::ExceptionBase & ExcNotActive()
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcMatrixDoesNotMatch()
#define AssertIndexRange(index, range)
#define DeclExceptionMsg(Exception, defaulttext)
static ::ExceptionBase & ExcInvalidObject()
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcCantCompareIterators()
static ::ExceptionBase & ExcVectorDoesNotMatch()
static ::ExceptionBase & ExcMessage(std::string arg1)
@ past_the_end
Iterator reached end of container.
types::fe_index get_fe_index_or_default(const DoFAccessor< structdim, dim, spacedim, level_dof_access > &cell, const types::fe_index fe_index)
void get_cell_dof_indices(const ::DoFCellAccessor< dim, spacedim, level_dof_access > &accessor, Implementation::dof_index_vector_type &dof_indices, const unsigned int fe_index)
TriaIterator< ::DoFAccessor< dim - 1, dim, spacedim, level_dof_access > > get_face(const ::DoFCellAccessor< dim, spacedim, level_dof_access > &cell, const unsigned int i, const std::integral_constant< int, 1 >)
constexpr types::global_dof_index invalid_dof_index
constexpr types::geometric_orientation default_geometric_orientation
constexpr types::fe_index invalid_fe_index
::VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
unsigned short int fe_index
void process_vertex_dofs(DoFHandler< dim, spacedim > &dof_handler, const unsigned int vertex_index, const types::fe_index fe_index, types::global_dof_index *&dof_indices_ptr, const DoFProcessor &dof_processor) const
void process_dofs(const DoFHandler< dim, spacedim > &dof_handler, const unsigned int obj_level, const unsigned int obj_index, const types::fe_index fe_index, const DoFMapping &mapping, const std::integral_constant< int, structdim >, types::global_dof_index *&dof_indices_ptr, const DoFProcessor &dof_processor) const
void process_dofs(const DoFHandler< dim, spacedim > &dof_handler, const unsigned int, const unsigned int obj_index, const types::fe_index fe_index, const DoFMapping &mapping, const std::integral_constant< int, structdim >, types::global_dof_index *&dof_indices_ptr, const DoFProcessor &dof_processor) const
MGDoFIndexProcessor(const unsigned int level)
void process_vertex_dofs(DoFHandler< dim, spacedim > &dof_handler, const unsigned int vertex_index, const types::fe_index, types::global_dof_index *&dof_indices_ptr, const DoFProcessor &dof_processor) const
static types::global_dof_index * get_array_ptr(const std::tuple<> &)
static void process_dof_indices(const ::DoFAccessor< structdim, dim, spacedim, level_dof_access > &accessor, const DoFIndicesType &const_dof_indices, const types::fe_index fe_index_, const DoFOperation &dof_operation, const DoFProcessor &dof_processor, const bool count_level_dofs)
static void get_dof_indices(const ::DoFAccessor< structdim, dim, spacedim, level_dof_access > &accessor, std::vector< types::global_dof_index > &dof_indices, const types::fe_index fe_index)
static void extract_subvector_to(const LinearAlgebra::TpetraWrappers::Vector< Number, MemorySpace > &values, const types::global_dof_index *cache_begin, const types::global_dof_index *cache_end, ForwardIterator local_values_begin)
static std::set< types::fe_index > get_active_fe_indices(const DoFHandler< dim, spacedim > &dof_handler, const unsigned int obj_level, const unsigned int obj_index, const std::integral_constant< int, structdim > &t)
static types::global_dof_index * get_array_ptr(const ArrayType &array)
static unsigned int get_array_length(const std::tuple<> &)
static unsigned int n_active_fe_indices(const DoFHandler< dim, spacedim > &dof_handler, const unsigned int obj_level, const unsigned int obj_index, const std::integral_constant< int, structdim > &)
static void set_dof_index(const DoFHandler< dim, spacedim > &dof_handler, const unsigned int obj_level, const unsigned int obj_index, const types::fe_index fe_index, const unsigned int local_index, const std::integral_constant< int, structdim > &dd, const types::global_dof_index global_index)
static void extract_subvector_to(const LinearAlgebra::EpetraWrappers::Vector &values, const types::global_dof_index *cache_begin, const types::global_dof_index *cache_end, ForwardIterator local_values_begin)
static void get_mg_dof_indices(const ::DoFAccessor< structdim, dim, spacedim, level_dof_access > &accessor, const int level, std::vector< types::global_dof_index > &dof_indices, const types::fe_index fe_index)
static bool fe_index_is_active(const DoFHandler< dim, spacedim > &dof_handler, const unsigned int obj_level, const unsigned int obj_index, const types::fe_index fe_index, const std::integral_constant< int, structdim > &)
static types::global_dof_index & get_mg_dof_index(const DoFHandler< dim, spacedim > &dof_handler, const std::unique_ptr< internal::DoFHandlerImplementation::DoFLevel< dim > > &, const std::unique_ptr< internal::DoFHandlerImplementation::DoFFaces< dim > > &mg_faces, const unsigned int obj_index, const types::fe_index fe_index, const unsigned int local_index, const std::integral_constant< int, 1 >)
static types::global_dof_index & get_mg_dof_index(const DoFHandler< dim, spacedim > &dof_handler, const std::unique_ptr< internal::DoFHandlerImplementation::DoFLevel< dim > > &mg_level, const std::unique_ptr< internal::DoFHandlerImplementation::DoFFaces< dim > > &, const unsigned int obj_index, const types::fe_index fe_index, const unsigned int local_index, const std::integral_constant< int, dim >)
static std::pair< unsigned int, unsigned int > process_object_range(const ::DoFAccessor< structdim, dim, spacedim, level_dof_access > accessor, const types::fe_index fe_index)
static types::global_dof_index & mg_vertex_dof_index(DoFHandler< dim, spacedim > &dof_handler, const int level, const unsigned int vertex_index, const unsigned int i)
static void set_mg_dof_indices(const ::DoFAccessor< structdim, dim, spacedim, level_dof_access > &accessor, const int level, const std::vector< types::global_dof_index > &dof_indices, const types::fe_index fe_index)
boost::container::small_vector<::types::global_dof_index, 27 > dof_index_vector_type
static void set_dof_indices(const ::DoFAccessor< structdim, dim, spacedim, level_dof_access > &accessor, const std::vector< types::global_dof_index > &dof_indices, const types::fe_index fe_index)
static void process_dof_index(const DoFHandler< dim, spacedim > &dof_handler, const unsigned int obj_level, const unsigned int obj_index, const types::fe_index fe_index, const unsigned int local_index, const std::integral_constant< int, structdim > &, GlobalIndexType &global_index, const DoFPProcessor &process)
static types::global_dof_index & get_mg_dof_index(const DoFHandler< 3, spacedim > &dof_handler, const std::unique_ptr< internal::DoFHandlerImplementation::DoFLevel< 3 > > &, const std::unique_ptr< internal::DoFHandlerImplementation::DoFFaces< 3 > > &mg_faces, const unsigned int obj_index, const types::fe_index fe_index, const unsigned int local_index, const std::integral_constant< int, 2 >)
static std::pair< unsigned int, unsigned int > process_object_range(const DoFHandler< dim, spacedim > &dof_handler, const unsigned int obj_level, const unsigned int obj_index, const types::fe_index fe_index, const std::integral_constant< int, structdim > &)
static void extract_subvector_to(const InputVector &values, const types::global_dof_index *cache, const types::global_dof_index *cache_end, ForwardIterator local_values_begin)
static std::vector< unsigned int > sort_indices(const types::global_dof_index *v_begin, const types::global_dof_index *v_end)
static types::fe_index nth_active_fe_index(const DoFHandler< dim, spacedim > &dof_handler, const unsigned int obj_level, const unsigned int obj_index, const unsigned int local_index, const std::integral_constant< int, structdim > &)
static void process_object(const DoFHandler< dim, spacedim > &dof_handler, const unsigned int obj_level, const unsigned int obj_index, const types::fe_index fe_index, const DoFMapping &mapping, const std::integral_constant< int, structdim > &dd, types::global_dof_index *&dof_indices_ptr, const DoFProcessor &process)
static types::global_dof_index get_dof_index(const DoFHandler< dim, spacedim > &dof_handler, const unsigned int obj_level, const unsigned int obj_index, const types::fe_index fe_index, const unsigned int local_index, const std::integral_constant< int, structdim > &dd)
static unsigned int n_dof_indices(const ::DoFAccessor< structdim, dim, spacedim, level_dof_access > &accessor, const types::fe_index fe_index_, const bool count_level_dofs)
static std::pair< unsigned int, unsigned int > process_object_range(::DoFInvalidAccessor< structdim, dim, spacedim >, const unsigned int)
static unsigned int get_array_length(const ArrayType &array)
static types::fe_index future_fe_index(const DoFCellAccessor< dim, spacedim, level_dof_access > &accessor)
static void clear_future_fe_index(const DoFCellAccessor< dim, spacedim, level_dof_access > &accessor)
static bool future_fe_index_set(const DoFCellAccessor< dim, spacedim, level_dof_access > &accessor)
static void set_active_fe_index(const DoFCellAccessor< dim, spacedim, level_dof_access > &accessor, const types::fe_index i)
static void set_future_fe_index(const DoFCellAccessor< dim, spacedim, level_dof_access > &accessor, const types::fe_index i)
static types::fe_index active_fe_index(const DoFCellAccessor< dim, spacedim, level_dof_access > &accessor)