Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
constraint_info.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2022 - 2023 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16
17#ifndef dealii_matrix_free_constraint_info_h
18#define dealii_matrix_free_constraint_info_h
19
20
21#include <deal.II/base/config.h>
22
24
29
30#include <limits>
31
33
34namespace internal
35{
36 namespace MatrixFreeFunctions
37 {
42 template <typename Number>
44 {
46
54 template <typename number2>
55 unsigned short
57 const std::vector<std::pair<types::global_dof_index, number2>>
58 &entries);
59
63 std::size_t
65
66 std::vector<std::pair<types::global_dof_index, double>>
68 std::vector<types::global_dof_index> constraint_indices;
69
70 std::pair<std::vector<Number>, types::global_dof_index> next_constraint;
71 std::map<std::vector<Number>,
75 };
76
77
78
84 template <int dim, typename Number>
86 {
87 public:
92 void
93 reinit(const DoFHandler<dim> &dof_handler,
94 const unsigned int n_cells,
95 const bool use_fast_hanging_node_algorithm = true);
96
97 void
99 const unsigned int cell_no,
100 const unsigned int mg_level,
102 const ::AffineConstraints<typename Number::value_type>
103 & constraints,
104 const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner);
105
109 void
110 reinit(const unsigned int n_cells);
111
112 void
114 const unsigned int cell_no,
115 const std::vector<types::global_dof_index> & dof_indices,
116 const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner);
117
118 void
120
121 template <typename T, typename VectorType>
122 void
123 read_write_operation(const T & operation,
124 VectorType & global_vector,
125 Number * local_vector,
126 const unsigned int first_cell,
127 const unsigned int n_cells,
128 const unsigned int n_dofs_per_cell,
129 const bool apply_constraints) const;
130
131 void
133 const unsigned int first_cell,
134 const unsigned int n_lanes_filled,
135 const bool transpose,
136 AlignedVector<Number> &evaluation_data_coarse) const;
137
141 std::size_t
143
144 private:
145 // for setup
147 std::vector<std::vector<unsigned int>> dof_indices_per_cell;
148 std::vector<std::vector<unsigned int>> plain_dof_indices_per_cell;
149 std::vector<std::vector<std::pair<unsigned short, unsigned short>>>
151
152 std::unique_ptr<HangingNodes<dim>> hanging_nodes;
153 std::vector<std::vector<unsigned int>> lexicographic_numbering;
154
155 public:
156 // for read_write_operation()
157 std::vector<unsigned int> dof_indices;
158 std::vector<std::pair<unsigned short, unsigned short>>
160 std::vector<std::pair<unsigned int, unsigned int>> row_starts;
161
162 std::vector<unsigned int> plain_dof_indices;
163 std::vector<unsigned int> row_starts_plain_indices;
164
165 // for constraint_pool_begin/end()
166 std::vector<typename Number::value_type> constraint_pool_data;
167 std::vector<unsigned int> constraint_pool_row_index;
168
169 std::vector<ShapeInfo<Number>> shape_infos;
170 std::vector<compressed_constraint_kind> hanging_node_constraint_masks;
171 std::vector<unsigned int> active_fe_indices;
172
173 private:
174 inline const typename Number::value_type *
175 constraint_pool_begin(const unsigned int row) const;
176
177 inline const typename Number::value_type *
178 constraint_pool_end(const unsigned int row) const;
179 };
180
181
182
183 // ------------------------- inline functions --------------------------
184
185 // NOLINTNEXTLINE(modernize-use-equals-default)
186 template <typename Number>
188 : constraints(FloatingPointComparator<Number>(
189 1. * std::numeric_limits<double>::epsilon() * 1024.))
190 {}
191
192
193
194 template <typename Number>
195 template <typename number2>
196 unsigned short
198 const std::vector<std::pair<types::global_dof_index, number2>> &entries)
199 {
200 next_constraint.first.resize(entries.size());
201 if (entries.size() > 0)
202 {
203 constraint_indices.resize(entries.size());
204 // Use assign so that values for nonmatching Number / number2 are
205 // converted:
206 constraint_entries.assign(entries.begin(), entries.end());
207 std::sort(constraint_entries.begin(),
208 constraint_entries.end(),
209 [](const std::pair<types::global_dof_index, double> &p1,
210 const std::pair<types::global_dof_index, double> &p2) {
211 return p1.second < p2.second;
212 });
213 for (types::global_dof_index j = 0; j < constraint_entries.size();
214 j++)
215 {
216 // copy the indices of the constraint entries after sorting.
217 constraint_indices[j] = constraint_entries[j].first;
218
219 // one_constraint takes the weights of the constraint
220 next_constraint.first[j] = constraint_entries[j].second;
221 }
222 }
223
224 // check whether or not constraint is already in pool. the initial
225 // implementation computed a hash value based on the truncated array (to
226 // given accuracy around 1e-13) in order to easily detect different
227 // arrays and then made a fine-grained check when the hash values were
228 // equal. this was quite lengthy and now we use a std::map with a
229 // user-defined comparator to compare floating point arrays to a
230 // tolerance 1e-13.
232 const auto position = constraints.find(next_constraint.first);
233 if (position != constraints.end())
234 insert_position = position->second;
235 else
236 {
237 next_constraint.second = constraints.size();
238 constraints.insert(next_constraint);
239 insert_position = next_constraint.second;
240 }
241
242 // we want to store the result as a short variable, so we have to make
243 // sure that the result does not exceed the limits when casting.
244 Assert(insert_position < (1 << (8 * sizeof(unsigned short))),
246 return static_cast<unsigned short>(insert_position);
247 }
248
249
250
251 template <int dim, typename Number>
252 inline void
254 const DoFHandler<dim> &dof_handler,
255 const unsigned int n_cells,
256 const bool use_fast_hanging_node_algorithm)
257 {
258 this->dof_indices_per_cell.resize(n_cells);
259 this->plain_dof_indices_per_cell.resize(n_cells);
260 this->constraint_indicator_per_cell.resize(n_cells);
261
262 // note: has_hanging_nodes() is a global operatrion
263 const bool has_hanging_nodes =
264 dof_handler.get_triangulation().has_hanging_nodes();
265
266 if (use_fast_hanging_node_algorithm && has_hanging_nodes)
267 {
268 hanging_nodes = std::make_unique<HangingNodes<dim>>(
269 dof_handler.get_triangulation());
270
271 hanging_node_constraint_masks.resize(n_cells);
272 }
273
274 const auto &fes = dof_handler.get_fe_collection();
275 lexicographic_numbering.resize(fes.size());
276 shape_infos.resize(fes.size());
277
278 for (unsigned int i = 0; i < fes.size(); ++i)
279 {
280 if (fes[i].reference_cell().is_hyper_cube())
281 {
282 const Quadrature<1> dummy_quadrature(
283 std::vector<Point<1>>(1, Point<1>()));
284 shape_infos[i].reinit(dummy_quadrature, fes[i], 0);
285 }
286 else
287 {
288 const auto dummy_quadrature =
289 fes[i].reference_cell().template get_gauss_type_quadrature<dim>(
290 1);
291 shape_infos[i].reinit(dummy_quadrature, fes[i], 0);
292 }
293
294 lexicographic_numbering[i] = shape_infos[i].lexicographic_numbering;
295 }
296 active_fe_indices.resize(n_cells);
297 }
298
299
300
301 template <int dim, typename Number>
302 inline void
303 ConstraintInfo<dim, Number>::reinit(const unsigned int n_cells)
304 {
305 this->dof_indices_per_cell.resize(n_cells);
306 this->plain_dof_indices_per_cell.resize(0);
307 this->constraint_indicator_per_cell.resize(n_cells);
308 }
309
310
311
312 template <int dim, typename Number>
313 inline void
315 const unsigned int cell_no,
316 const unsigned int mg_level,
318 const ::AffineConstraints<typename Number::value_type> &constraints,
319 const std::shared_ptr<const Utilities::MPI::Partitioner> & partitioner)
320 {
321 std::vector<types::global_dof_index> local_dof_indices(
322 cell->get_fe().n_dofs_per_cell());
323 std::vector<types::global_dof_index> local_dof_indices_lex(
324 cell->get_fe().n_dofs_per_cell());
325
326 if (mg_level == numbers::invalid_unsigned_int)
327 cell->get_dof_indices(local_dof_indices);
328 else
329 cell->get_mg_dof_indices(local_dof_indices);
330
331 {
332 AssertIndexRange(cell->active_fe_index(), shape_infos.size());
333
334 const auto &lexicographic_numbering =
335 shape_infos[cell->active_fe_index()].lexicographic_numbering;
336
337 AssertDimension(lexicographic_numbering.size(),
338 local_dof_indices.size());
339
340 for (unsigned int i = 0; i < cell->get_fe().n_dofs_per_cell(); ++i)
341 local_dof_indices_lex[i] =
342 local_dof_indices[lexicographic_numbering[i]];
343 }
344
345 std::pair<unsigned short, unsigned short> constraint_iterator(0, 0);
346
347 AssertIndexRange(cell_no, this->constraint_indicator_per_cell.size());
348 AssertIndexRange(cell_no, this->dof_indices_per_cell.size());
349 AssertIndexRange(cell_no, this->plain_dof_indices_per_cell.size());
350 AssertIndexRange(cell_no, this->plain_dof_indices_per_cell.size());
351
352 auto &constraint_indicator = this->constraint_indicator_per_cell[cell_no];
353 auto &dof_indices = this->dof_indices_per_cell[cell_no];
354 auto &plain_dof_indices = this->plain_dof_indices_per_cell[cell_no];
355
356 AssertDimension(constraint_indicator_per_cell[cell_no].size(), 0);
357 AssertDimension(dof_indices_per_cell[cell_no].size(), 0);
358 AssertDimension(plain_dof_indices_per_cell[cell_no].size(), 0);
359
360 const auto global_to_local =
361 [&](const types::global_dof_index global_index) -> unsigned int {
362 if (partitioner)
363 return partitioner->global_to_local(global_index);
364 else
365 return global_index;
366 };
367
368 // plain indices
369 plain_dof_indices.resize(local_dof_indices_lex.size());
370 for (unsigned int i = 0; i < local_dof_indices_lex.size(); ++i)
371 plain_dof_indices[i] = global_to_local(local_dof_indices_lex[i]);
372
373 if (hanging_nodes)
374 {
375 AssertIndexRange(cell_no, this->hanging_node_constraint_masks.size());
376 AssertIndexRange(cell_no, this->active_fe_indices.size());
377
378 std::vector<ConstraintKinds> mask(cell->get_fe().n_components());
379 hanging_nodes->setup_constraints(
380 cell, {}, lexicographic_numbering, local_dof_indices_lex, mask);
381
382 hanging_node_constraint_masks[cell_no] = compress(mask[0], dim);
383 active_fe_indices[cell_no] = cell->active_fe_index();
384 }
385
386 for (auto current_dof : local_dof_indices_lex)
387 {
388 const auto *entries_ptr =
389 constraints.get_constraint_entries(current_dof);
390
391 // dof is constrained
392 if (entries_ptr != nullptr)
393 {
394 const auto & entries = *entries_ptr;
395 const types::global_dof_index n_entries = entries.size();
396 if (n_entries == 1 &&
397 std::abs(entries[0].second - 1.) <
398 100 * std::numeric_limits<double>::epsilon())
399 {
400 current_dof = entries[0].first;
401 goto no_constraint;
402 }
403
404 constraint_indicator.push_back(constraint_iterator);
405 constraint_indicator.back().second =
406 constraint_values.insert_entries(entries);
407
408 // reset constraint iterator for next round
409 constraint_iterator.first = 0;
410
411 if (n_entries > 0)
412 {
413 const std::vector<types::global_dof_index>
414 &constraint_indices = constraint_values.constraint_indices;
415 for (unsigned int j = 0; j < n_entries; ++j)
416 {
417 dof_indices.push_back(
418 global_to_local(constraint_indices[j]));
419 }
420 }
421 }
422 else
423 {
424 no_constraint:
425 dof_indices.push_back(global_to_local(current_dof));
426
427 // make sure constraint_iterator.first is always within the
428 // bounds of unsigned short
429 Assert(constraint_iterator.first <
430 (1 << (8 * sizeof(unsigned short))) - 1,
432 constraint_iterator.first++;
433 }
434 }
435 }
436
437
438
439 template <int dim, typename Number>
440 inline void
442 const unsigned int cell_no,
443 const std::vector<types::global_dof_index> &local_dof_indices_lex,
444 const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner)
445 {
446 const auto global_to_local =
447 [&](const types::global_dof_index global_index) -> unsigned int {
448 if (partitioner)
449 return partitioner->global_to_local(global_index);
450 else
451 return global_index;
452 };
453
454 std::pair<unsigned short, unsigned short> constraint_iterator(0, 0);
455
456 auto &constraint_indicator = this->constraint_indicator_per_cell[cell_no];
457 auto &dof_indices = this->dof_indices_per_cell[cell_no];
458
459 for (const auto current_dof : local_dof_indices_lex)
460 {
461 // dof is constrained
462 if (current_dof == numbers::invalid_dof_index)
463 {
464 const std::vector<
465 std::pair<types::global_dof_index, typename Number::value_type>>
466 entries = {};
467
468 constraint_indicator.push_back(constraint_iterator);
469 constraint_indicator.back().second =
470 constraint_values.insert_entries(entries);
471
472 constraint_iterator.first = 0;
473 }
474 else
475 {
476 dof_indices.push_back(global_to_local(current_dof));
477
478 // make sure constraint_iterator.first is always within the
479 // bounds of unsigned short
480 Assert(constraint_iterator.first <
481 (1 << (8 * sizeof(unsigned short))) - 1,
483 constraint_iterator.first++;
484 }
485 }
486 }
487
488
489
490 template <int dim, typename Number>
491 inline void
493 {
494 this->dof_indices = {};
495 this->plain_dof_indices = {};
496 this->constraint_indicator = {};
497
498 this->row_starts = {};
499 this->row_starts.emplace_back(0, 0);
500
501 if (plain_dof_indices_per_cell.empty() == false)
502 {
503 this->row_starts_plain_indices = {};
504 this->row_starts_plain_indices.emplace_back(0);
505 }
506
507 for (unsigned int i = 0; i < dof_indices_per_cell.size(); ++i)
508 {
509 this->dof_indices.insert(this->dof_indices.end(),
510 dof_indices_per_cell[i].begin(),
511 dof_indices_per_cell[i].end());
512 this->constraint_indicator.insert(
513 this->constraint_indicator.end(),
514 constraint_indicator_per_cell[i].begin(),
515 constraint_indicator_per_cell[i].end());
516
517 this->row_starts.emplace_back(this->dof_indices.size(),
518 this->constraint_indicator.size());
519
520 if (plain_dof_indices_per_cell.empty() == false)
521 {
522 this->plain_dof_indices.insert(
523 this->plain_dof_indices.end(),
524 plain_dof_indices_per_cell[i].begin(),
525 plain_dof_indices_per_cell[i].end());
526
527 this->row_starts_plain_indices.emplace_back(
528 this->plain_dof_indices.size());
529 }
530 }
531
532 std::vector<const std::vector<double> *> constraints(
533 constraint_values.constraints.size());
534 unsigned int length = 0;
535 for (const auto &it : constraint_values.constraints)
536 {
537 AssertIndexRange(it.second, constraints.size());
538 constraints[it.second] = &it.first;
539 length += it.first.size();
540 }
541
542 constraint_pool_data.clear();
543 constraint_pool_data.reserve(length);
544 constraint_pool_row_index.reserve(constraint_values.constraints.size() +
545 1);
546 constraint_pool_row_index.resize(1, 0);
547
548 for (const auto &constraint : constraints)
549 {
550 Assert(constraint != nullptr, ExcInternalError());
551 constraint_pool_data.insert(constraint_pool_data.end(),
552 constraint->begin(),
553 constraint->end());
554 constraint_pool_row_index.push_back(constraint_pool_data.size());
555 }
556
557 AssertDimension(constraint_pool_data.size(), length);
558
559 dof_indices_per_cell.clear();
560 constraint_indicator_per_cell.clear();
561
562 if (hanging_nodes &&
563 std::all_of(hanging_node_constraint_masks.begin(),
564 hanging_node_constraint_masks.end(),
565 [](const auto i) {
566 return i == unconstrained_compressed_constraint_kind;
567 }))
568 hanging_node_constraint_masks.clear();
569 }
570
571
572
573 template <int dim, typename Number>
574 template <typename T, typename VectorType>
575 inline void
577 const T & operation,
578 VectorType & global_vector,
579 Number * local_vector,
580 const unsigned int first_cell,
581 const unsigned int n_cells,
582 const unsigned int n_dofs_per_cell,
583 const bool apply_constraints) const
584 {
585 if ((row_starts_plain_indices.empty() == false) &&
586 (apply_constraints == false))
587 {
588 for (unsigned int v = 0; v < n_cells; ++v)
589 {
590 const unsigned int cell_index = first_cell + v;
591 const unsigned int *dof_indices =
592 this->plain_dof_indices.data() +
593 this->row_starts_plain_indices[cell_index];
594
595 for (unsigned int i = 0; i < n_dofs_per_cell; ++dof_indices, ++i)
596 operation.process_dof(*dof_indices,
597 global_vector,
598 local_vector[i][v]);
599 }
600
601 return;
602 }
603
604 for (unsigned int v = 0; v < n_cells; ++v)
605 {
606 const unsigned int cell_index = first_cell + v;
607 const unsigned int *dof_indices =
608 this->dof_indices.data() + this->row_starts[cell_index].first;
609 unsigned int index_indicators = this->row_starts[cell_index].second;
610 unsigned int next_index_indicators =
611 this->row_starts[cell_index + 1].second;
612
613 unsigned int ind_local = 0;
614 for (; index_indicators != next_index_indicators; ++index_indicators)
615 {
616 const std::pair<unsigned short, unsigned short> indicator =
617 this->constraint_indicator[index_indicators];
618
619 // run through values up to next constraint
620 for (unsigned int j = 0; j < indicator.first; ++j)
621 operation.process_dof(dof_indices[j],
622 global_vector,
623 local_vector[ind_local + j][v]);
624
625 ind_local += indicator.first;
626 dof_indices += indicator.first;
627
628 // constrained case: build the local value as a linear
629 // combination of the global value according to constraints
630 typename Number::value_type value;
631 operation.pre_constraints(local_vector[ind_local][v], value);
632
633 const typename Number::value_type *data_val =
634 this->constraint_pool_begin(indicator.second);
635 const typename Number::value_type *end_pool =
636 this->constraint_pool_end(indicator.second);
637 for (; data_val != end_pool; ++data_val, ++dof_indices)
638 operation.process_constraint(*dof_indices,
639 *data_val,
640 global_vector,
641 value);
642
643 operation.post_constraints(value, local_vector[ind_local][v]);
644 ind_local++;
645 }
646
647 AssertIndexRange(ind_local, n_dofs_per_cell + 1);
648
649 for (; ind_local < n_dofs_per_cell; ++dof_indices, ++ind_local)
650 operation.process_dof(*dof_indices,
651 global_vector,
652 local_vector[ind_local][v]);
653 }
654 }
655
656
657
658 template <int dim, typename Number>
659 inline void
661 const unsigned int first_cell,
662 const unsigned int n_lanes_filled,
663 const bool transpose,
664 AlignedVector<Number> &evaluation_data_coarse) const
665 {
666 if (hanging_node_constraint_masks.size() == 0)
667 return;
668
670 Number::size()>
671 constraint_mask;
672
673 bool hn_available = false;
674
675 for (unsigned int v = 0; v < n_lanes_filled; ++v)
676 {
677 const auto mask = hanging_node_constraint_masks[first_cell + v];
678
679 constraint_mask[v] = mask;
680
682 }
683
684 if (hn_available == true)
685 {
686 for (unsigned int v = n_lanes_filled; v < Number::size(); ++v)
687 constraint_mask[v] = unconstrained_compressed_constraint_kind;
688
689 for (unsigned int i = 1; i < n_lanes_filled; ++i)
690 AssertDimension(active_fe_indices[first_cell],
691 active_fe_indices[first_cell + i]);
692
693 const auto &shape_info = shape_infos[active_fe_indices[first_cell]];
694
696 dim,
697 typename Number::value_type,
698 Number>::apply(shape_info.n_components,
699 shape_info.data.front().fe_degree,
700 shape_info,
701 transpose,
702 constraint_mask,
703 evaluation_data_coarse.begin());
704 }
705 }
706
707
708
709 template <int dim, typename Number>
710 inline const typename Number::value_type *
712 const unsigned int row) const
713 {
714 AssertIndexRange(row, constraint_pool_row_index.size() - 1);
715 return constraint_pool_data.empty() ?
716 nullptr :
717 constraint_pool_data.data() + constraint_pool_row_index[row];
718 }
719
720
721
722 template <int dim, typename Number>
723 inline const typename Number::value_type *
725 const unsigned int row) const
726 {
727 AssertIndexRange(row, constraint_pool_row_index.size() - 1);
728 return constraint_pool_data.empty() ?
729 nullptr :
730 constraint_pool_data.data() + constraint_pool_row_index[row + 1];
731 }
732
733
734
735 template <int dim, typename Number>
736 inline std::size_t
738 {
739 std::size_t size = 0;
740
741 size += MemoryConsumption::memory_consumption(constraint_values);
742 size += MemoryConsumption::memory_consumption(dof_indices_per_cell);
743 size += MemoryConsumption::memory_consumption(plain_dof_indices_per_cell);
744 size +=
745 MemoryConsumption::memory_consumption(constraint_indicator_per_cell);
746
747 if (hanging_nodes)
748 size += MemoryConsumption::memory_consumption(*hanging_nodes);
749
750 size += MemoryConsumption::memory_consumption(lexicographic_numbering);
751 size += MemoryConsumption::memory_consumption(dof_indices);
752 size += MemoryConsumption::memory_consumption(constraint_indicator);
753 size += MemoryConsumption::memory_consumption(row_starts);
754 size += MemoryConsumption::memory_consumption(plain_dof_indices);
755 size += MemoryConsumption::memory_consumption(row_starts_plain_indices);
756 size += MemoryConsumption::memory_consumption(constraint_pool_data);
757 size += MemoryConsumption::memory_consumption(constraint_pool_row_index);
758 size += MemoryConsumption::memory_consumption(shape_infos);
759 size +=
760 MemoryConsumption::memory_consumption(hanging_node_constraint_masks);
761 size += MemoryConsumption::memory_consumption(active_fe_indices);
762
763 return size;
764 }
765
766
767
768 template <typename Number>
769 inline std::size_t
771 {
772 std::size_t size = 0;
773
774 size += MemoryConsumption::memory_consumption(constraint_entries);
775 size += MemoryConsumption::memory_consumption(constraint_indices);
776
777 // TODO: map does not have memory_consumption()
778 // size += MemoryConsumption::memory_consumption(constraints);
779
780 return size;
781 }
782
783 } // namespace MatrixFreeFunctions
784} // namespace internal
785
787
788#endif
iterator begin()
const hp::FECollection< dim, spacedim > & get_fe_collection() const
const Triangulation< dim, spacedim > & get_triangulation() const
Definition point.h:112
virtual bool has_hanging_nodes() const
void read_dof_indices(const unsigned int cell_no, const std::vector< types::global_dof_index > &dof_indices, const std::shared_ptr< const Utilities::MPI::Partitioner > &partitioner)
std::vector< std::vector< unsigned int > > lexicographic_numbering
void read_write_operation(const T &operation, VectorType &global_vector, Number *local_vector, const unsigned int first_cell, const unsigned int n_cells, const unsigned int n_dofs_per_cell, const bool apply_constraints) const
std::vector< std::pair< unsigned short, unsigned short > > constraint_indicator
void reinit(const DoFHandler< dim > &dof_handler, const unsigned int n_cells, const bool use_fast_hanging_node_algorithm=true)
std::vector< std::vector< std::pair< unsigned short, unsigned short > > > constraint_indicator_per_cell
std::vector< unsigned int > constraint_pool_row_index
const Number::value_type * constraint_pool_begin(const unsigned int row) const
std::vector< std::vector< unsigned int > > plain_dof_indices_per_cell
std::vector< std::vector< unsigned int > > dof_indices_per_cell
std::vector< ShapeInfo< Number > > shape_infos
std::unique_ptr< HangingNodes< dim > > hanging_nodes
std::vector< unsigned int > row_starts_plain_indices
const Number::value_type * constraint_pool_end(const unsigned int row) const
std::vector< compressed_constraint_kind > hanging_node_constraint_masks
void apply_hanging_node_constraints(const unsigned int first_cell, const unsigned int n_lanes_filled, const bool transpose, AlignedVector< Number > &evaluation_data_coarse) const
void read_dof_indices(const unsigned int cell_no, const unsigned int mg_level, const TriaIterator< DoFCellAccessor< dim, dim, false > > &cell, const ::AffineConstraints< typename Number::value_type > &constraints, const std::shared_ptr< const Utilities::MPI::Partitioner > &partitioner)
std::vector< std::pair< unsigned int, unsigned int > > row_starts
std::vector< typename Number::value_type > constraint_pool_data
void reinit(const unsigned int n_cells)
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
DerivativeForm< 1, spacedim, dim, Number > transpose(const DerivativeForm< 1, dim, spacedim, Number > &DF)
Point< 2 > second
Definition grid_out.cc:4616
unsigned int cell_index
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
std::enable_if_t< std::is_fundamental< T >::value, std::size_t > memory_consumption(const T &t)
compressed_constraint_kind compress(const ConstraintKinds kind_in, const unsigned int dim)
std::uint8_t compressed_constraint_kind
Definition dof_info.h:83
constexpr compressed_constraint_kind unconstrained_compressed_constraint_kind
static const unsigned int invalid_unsigned_int
Definition types.h:213
const types::global_dof_index invalid_dof_index
Definition types.h:233
STL namespace.
::VectorizedArray< Number, width > abs(const ::VectorizedArray< Number, width > &)
unsigned int global_dof_index
Definition types.h:82
std::vector< types::global_dof_index > constraint_indices
std::map< std::vector< Number >, types::global_dof_index, FloatingPointComparator< Number > > constraints
unsigned short insert_entries(const std::vector< std::pair< types::global_dof_index, number2 > > &entries)
std::vector< std::pair< types::global_dof_index, double > > constraint_entries
std::pair< std::vector< Number >, types::global_dof_index > next_constraint