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
mapping_q_cache.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2019 - 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
19
21
22#include <deal.II/fe/fe_dgq.h>
24#include <deal.II/fe/fe_q.h>
25#include <deal.II/fe/fe_tools.h>
28
33
34#include <functional>
35
37
38template <int dim, int spacedim>
40 const unsigned int polynomial_degree)
41 : MappingQ<dim, spacedim>(polynomial_degree)
42 , uses_level_info(false)
43{}
44
45
46
47template <int dim, int spacedim>
49 const MappingQCache<dim, spacedim> &mapping)
50 : MappingQ<dim, spacedim>(mapping)
51 , support_point_cache(mapping.support_point_cache)
52 , uses_level_info(mapping.uses_level_info)
53{}
54
55
56
57template <int dim, int spacedim>
59{
60 // When this object goes out of scope, we want the cache to get cleared and
61 // free its memory before the signal is disconnected in order to not work on
62 // invalid memory that has been left back by freeing an object of this
63 // class.
64 support_point_cache.reset();
65 clear_signal.disconnect();
66}
67
68
69
70template <int dim, int spacedim>
71std::unique_ptr<Mapping<dim, spacedim>>
73{
74 return std::make_unique<MappingQCache<dim, spacedim>>(*this);
75}
76
77
78
79template <int dim, int spacedim>
80bool
82{
83 return false;
84}
85
86
87
88template <int dim, int spacedim>
89void
91 const Mapping<dim, spacedim> & mapping,
93{
94 // FE and FEValues in the case they are needed
97 fe_values_all;
98
99 this->initialize(
101 [&](const typename Triangulation<dim, spacedim>::cell_iterator &cell) {
102 const auto mapping_q =
103 dynamic_cast<const MappingQ<dim, spacedim> *>(&mapping);
104 if (mapping_q != nullptr && this->get_degree() == mapping_q->get_degree())
105 {
106 return mapping_q->compute_mapping_support_points(cell);
107 }
108 else
109 {
110 // get FEValues (thread-safe); in the case that this thread has not
111 // created a an FEValues object yet, this helper-function also
112 // creates one with the right quadrature rule
113 auto &fe_values = fe_values_all.get();
114 if (fe_values.get() == nullptr)
115 {
116 QGaussLobatto<dim> quadrature_gl(this->polynomial_degree + 1);
117
118 std::vector<Point<dim>> quadrature_points;
119 for (const auto i :
120 FETools::hierarchic_to_lexicographic_numbering<dim>(
121 this->polynomial_degree))
122 quadrature_points.push_back(quadrature_gl.point(i));
123 Quadrature<dim> quadrature(quadrature_points);
124
125 fe_values = std::make_unique<FEValues<dim, spacedim>>(
126 mapping, fe, quadrature, update_quadrature_points);
127 }
128
129 fe_values->reinit(cell);
130 return fe_values->get_quadrature_points();
131 }
132 });
133}
134
135
136
137template <int dim, int spacedim>
138void
141 const MappingQ<dim, spacedim> & mapping)
142{
143 this->initialize(mapping, triangulation);
144}
145
146
147
148template <int dim, int spacedim>
149void
152 const std::function<std::vector<Point<spacedim>>(
154 &compute_points_on_cell)
155{
156 clear_signal.disconnect();
157 clear_signal = triangulation.signals.any_change.connect(
158 [&]() -> void { this->support_point_cache.reset(); });
159
160 support_point_cache =
161 std::make_shared<std::vector<std::vector<std::vector<Point<spacedim>>>>>(
162 triangulation.n_levels());
163 for (unsigned int l = 0; l < triangulation.n_levels(); ++l)
164 (*support_point_cache)[l].resize(triangulation.n_raw_cells(l));
165
167 triangulation.begin(),
168 triangulation.end(),
169 [&](const typename Triangulation<dim, spacedim>::cell_iterator &cell,
170 void *,
171 void *) {
172 (*support_point_cache)[cell->level()][cell->index()] =
173 compute_points_on_cell(cell);
174 AssertDimension(
175 (*support_point_cache)[cell->level()][cell->index()].size(),
176 Utilities::pow(this->get_degree() + 1, dim));
177 },
178 /* copier */ std::function<void(void *)>(),
179 /* scratch_data */ nullptr,
180 /* copy_data */ nullptr,
182 /* chunk_size = */ 1);
183
184 uses_level_info = true;
185}
186
187
188
189template <int dim, int spacedim>
190void
192 const Mapping<dim, spacedim> & mapping,
194 const std::function<Point<spacedim>(
196 const Point<spacedim> &)> & transformation_function,
197 const bool function_describes_relative_displacement)
198{
199 // FE and FEValues in the case they are needed
202 fe_values_all;
203
204 this->initialize(
205 tria,
206 [&](const typename Triangulation<dim, spacedim>::cell_iterator &cell) {
207 std::vector<Point<spacedim>> points;
208
209 const auto mapping_q =
210 dynamic_cast<const MappingQ<dim, spacedim> *>(&mapping);
211
212 if (mapping_q != nullptr && this->get_degree() == mapping_q->get_degree())
213 {
214 points = mapping_q->compute_mapping_support_points(cell);
215 }
216 else
217 {
218 // get FEValues (thread-safe); in the case that this thread has not
219 // created a an FEValues object yet, this helper-function also
220 // creates one with the right quadrature rule
221 auto &fe_values = fe_values_all.get();
222 if (fe_values.get() == nullptr)
223 {
224 QGaussLobatto<dim> quadrature_gl(this->polynomial_degree + 1);
225
226 std::vector<Point<dim>> quadrature_points;
227 for (const auto i :
228 FETools::hierarchic_to_lexicographic_numbering<dim>(
229 this->polynomial_degree))
230 quadrature_points.push_back(quadrature_gl.point(i));
231 Quadrature<dim> quadrature(quadrature_points);
232
233 fe_values = std::make_unique<FEValues<dim, spacedim>>(
234 mapping, fe, quadrature, update_quadrature_points);
235 }
236
237 fe_values->reinit(cell);
238 points = fe_values->get_quadrature_points();
239 }
240
241 for (auto &p : points)
242 if (function_describes_relative_displacement)
243 p += transformation_function(cell, p);
244 else
245 p = transformation_function(cell, p);
246
247 return points;
248 });
249
250 uses_level_info = true;
251}
252
253
254
255template <int dim, int spacedim>
256void
258 const Mapping<dim, spacedim> & mapping,
260 const Function<spacedim> & transformation_function,
261 const bool function_describes_relative_displacement)
262{
263 AssertDimension(transformation_function.n_components, spacedim);
264
265 this->initialize(
266 mapping,
267 tria,
268 [&](const auto &, const auto &point) {
269 Point<spacedim> new_point;
270 for (unsigned int c = 0; c < spacedim; ++c)
271 new_point[c] = transformation_function.value(point, c);
272 return new_point;
273 },
274 function_describes_relative_displacement);
275
276 uses_level_info = true;
277}
278
279
280
281namespace
282{
283 template <typename VectorType>
284 void
285 copy_locally_owned_data_from(
286 const VectorType &vector,
288 &vector_ghosted)
289 {
291 temp.reinit(vector.locally_owned_elements());
293 vector_ghosted.import_elements(temp, VectorOperation::insert);
294 }
295} // namespace
296
297
298
299template <int dim, int spacedim>
300template <typename VectorType>
301void
303 const Mapping<dim, spacedim> & mapping,
304 const DoFHandler<dim, spacedim> &dof_handler,
305 const VectorType & vector,
306 const bool vector_describes_relative_displacement)
307{
308 AssertDimension(dof_handler.get_fe_collection().size(), 1);
309 const FiniteElement<dim, spacedim> &fe = dof_handler.get_fe();
311 AssertDimension(fe.element_multiplicity(0), spacedim);
312
313 const unsigned int is_fe_q =
314 dynamic_cast<const FE_Q<dim, spacedim> *>(&fe.base_element(0)) != nullptr;
315 const unsigned int is_fe_dgq =
316 dynamic_cast<const FE_DGQ<dim, spacedim> *>(&fe.base_element(0)) != nullptr;
317
318 const auto lexicographic_to_hierarchic_numbering =
320 FETools::hierarchic_to_lexicographic_numbering<spacedim>(
321 this->get_degree()));
322
323 // Step 1: copy global vector so that the ghost values are such that the
324 // cache can be set up for all ghost cells
326 vector_ghosted;
327 IndexSet locally_relevant_dofs;
328 DoFTools::extract_locally_relevant_dofs(dof_handler, locally_relevant_dofs);
329 vector_ghosted.reinit(dof_handler.locally_owned_dofs(),
330 locally_relevant_dofs,
331 dof_handler.get_communicator());
332 copy_locally_owned_data_from(vector, vector_ghosted);
333 vector_ghosted.update_ghost_values();
334
335 // FE and FEValues in the case they are needed
336 FE_Nothing<dim, spacedim> fe_nothing;
338 fe_values_all;
339
340 // Interpolation of values is needed if we cannot just read off locations
341 // from the solution vectors (as in the case of FE_Q and FE_DGQ with the
342 // same polynomial degree as this class has).
343 const bool interpolation_of_values_is_needed =
344 ((is_fe_q || is_fe_dgq) && fe.degree == this->get_degree()) == false;
345
346 // Step 2: loop over all cells
347 this->initialize(
348 dof_handler.get_triangulation(),
349 [&](const typename Triangulation<dim, spacedim>::cell_iterator &cell_tria)
350 -> std::vector<Point<spacedim>> {
351 const bool is_active_non_artificial_cell =
352 (cell_tria->is_active() == true) &&
353 (cell_tria->is_artificial() == false);
354
355 const typename DoFHandler<dim, spacedim>::cell_iterator cell_dofs(
356 &cell_tria->get_triangulation(),
357 cell_tria->level(),
358 cell_tria->index(),
359 &dof_handler);
360
361 const auto mapping_q =
362 dynamic_cast<const MappingQ<dim, spacedim> *>(&mapping);
363
364 // Step 2a) set up and reinit FEValues (if needed)
365 if (
366 ((vector_describes_relative_displacement ||
367 (is_active_non_artificial_cell == false)) &&
368 ((mapping_q != nullptr &&
369 this->get_degree() == mapping_q->get_degree()) ==
370 false)) /*condition 1: points need to be computed via FEValues*/
371 ||
372 (is_active_non_artificial_cell && interpolation_of_values_is_needed) /*condition 2: interpolation of values is needed*/)
373 {
374 // get FEValues (thread-safe); in the case that this thread has
375 // not created a an FEValues object yet, this helper-function also
376 // creates one with the right quadrature rule
377 auto &fe_values = fe_values_all.get();
378 if (fe_values.get() == nullptr)
379 {
380 QGaussLobatto<dim> quadrature_gl(this->polynomial_degree + 1);
381
382 std::vector<Point<dim>> quadrature_points;
383 for (const auto i :
384 FETools::hierarchic_to_lexicographic_numbering<dim>(
385 this->polynomial_degree))
386 quadrature_points.push_back(quadrature_gl.point(i));
387 Quadrature<dim> quadrature(quadrature_points);
388
389 fe_values = std::make_unique<FEValues<dim, spacedim>>(
390 mapping,
391 interpolation_of_values_is_needed ?
392 fe :
393 static_cast<const FiniteElement<dim, spacedim> &>(fe_nothing),
394 quadrature,
395 update_quadrature_points | update_values);
396 }
397
398 if (interpolation_of_values_is_needed)
399 fe_values->reinit(cell_dofs);
400 else
401 fe_values->reinit(cell_tria);
402 }
403
404 std::vector<Point<spacedim>> result;
405
406 // Step 2b) read of quadrature points in the relative displacement case
407 // note: we also take this path for non-active or artificial cells so that
408 // these cells are filled with some useful data
409 if (vector_describes_relative_displacement ||
410 is_active_non_artificial_cell == false)
411 {
412 if (mapping_q != nullptr &&
413 this->get_degree() == mapping_q->get_degree())
414 result = mapping_q->compute_mapping_support_points(cell_tria);
415 else
416 result = fe_values_all.get()->get_quadrature_points();
417
418 // for non-active or artificial cells we are done here and return
419 // the absolute positions, since the provided vector cannot contain
420 // any useful information for these cells
421 if (is_active_non_artificial_cell == false)
422 return result;
423 }
424 else
425 {
426 result.resize(
427 Utilities::pow<unsigned int>(this->get_degree() + 1, dim));
428 }
429
430 // Step 2c) read global vector and adjust points accordingly
431 if (interpolation_of_values_is_needed == false)
432 {
433 // case 1: FE_Q or FE_DGQ with same degree as this class has; this
434 // is the simple case since no interpolation is needed
435 std::vector<types::global_dof_index> dof_indices(
436 fe.n_dofs_per_cell());
437 cell_dofs->get_dof_indices(dof_indices);
438
439 for (unsigned int i = 0; i < dof_indices.size(); ++i)
440 {
441 const auto id = fe.system_to_component_index(i);
442
443 if (is_fe_q)
444 {
445 // case 1a: FE_Q
446 if (vector_describes_relative_displacement)
447 result[id.second][id.first] +=
448 vector_ghosted(dof_indices[i]);
449 else
450 result[id.second][id.first] =
451 vector_ghosted(dof_indices[i]);
452 }
453 else
454 {
455 // case 1b: FE_DGQ
456 if (vector_describes_relative_displacement)
457 result[lexicographic_to_hierarchic_numbering[id.second]]
458 [id.first] += vector_ghosted(dof_indices[i]);
459 else
460 result[lexicographic_to_hierarchic_numbering[id.second]]
461 [id.first] = vector_ghosted(dof_indices[i]);
462 }
463 }
464 }
465 else
466 {
467 // case 2: general case; interpolation is needed
468 // note: the following code could be optimized for tensor-product
469 // elements via application of sum factorization as is done on
470 // MatrixFree/FEEvaluation
471 auto &fe_values = fe_values_all.get();
472
473 std::vector<Vector<typename VectorType::value_type>> values(
474 fe_values->n_quadrature_points,
476
477 fe_values->get_function_values(vector_ghosted, values);
478
479 for (unsigned int q = 0; q < fe_values->n_quadrature_points; ++q)
480 for (unsigned int c = 0; c < spacedim; ++c)
481 if (vector_describes_relative_displacement)
482 result[q][c] += values[q][c];
483 else
484 result[q][c] = values[q][c];
485 }
486
487 return result;
488 });
489
490 uses_level_info = false;
491}
492
493
494
495template <int dim, int spacedim>
496template <typename VectorType>
497void
499 const Mapping<dim, spacedim> & mapping,
500 const DoFHandler<dim, spacedim> &dof_handler,
501 const MGLevelObject<VectorType> &vectors,
502 const bool vector_describes_relative_displacement)
503{
504 AssertDimension(dof_handler.get_fe_collection().size(), 1);
505 const FiniteElement<dim, spacedim> &fe = dof_handler.get_fe();
507 AssertDimension(fe.element_multiplicity(0), spacedim);
508 AssertDimension(0, vectors.min_level());
510 vectors.max_level());
511
512 const unsigned int is_fe_q =
513 dynamic_cast<const FE_Q<dim, spacedim> *>(&fe.base_element(0)) != nullptr;
514 const unsigned int is_fe_dgq =
515 dynamic_cast<const FE_DGQ<dim, spacedim> *>(&fe.base_element(0)) != nullptr;
516
517 const auto lexicographic_to_hierarchic_numbering =
519 FETools::hierarchic_to_lexicographic_numbering<spacedim>(
520 this->get_degree()));
521
522 // Step 1: copy global vector so that the ghost values are such that the
523 // cache can be set up for all ghost cells
526 vectors_ghosted(vectors.min_level(), vectors.max_level());
527
528 for (unsigned int l = vectors.min_level(); l <= vectors.max_level(); ++l)
529 {
530 IndexSet locally_relevant_dofs;
532 l,
533 locally_relevant_dofs);
534 vectors_ghosted[l].reinit(dof_handler.locally_owned_mg_dofs(l),
535 locally_relevant_dofs,
536 dof_handler.get_communicator());
537 copy_locally_owned_data_from(vectors[l], vectors_ghosted[l]);
538 vectors_ghosted[l].update_ghost_values();
539 }
540
541 // FE and FEValues in the case they are needed
542 FE_Nothing<dim, spacedim> fe_nothing;
544 fe_values_all;
545
546 // Interpolation of values is needed if we cannot just read off locations
547 // from the solution vectors (as in the case of FE_Q and FE_DGQ with the
548 // same polynomial degree as this class has).
549 const bool interpolation_of_values_is_needed =
550 ((is_fe_q || is_fe_dgq) && fe.degree == this->get_degree()) == false;
551
552 // Step 2: loop over all cells
553 this->initialize(
554 dof_handler.get_triangulation(),
555 [&](const typename Triangulation<dim, spacedim>::cell_iterator &cell_tria)
556 -> std::vector<Point<spacedim>> {
557 const bool is_non_artificial_cell =
558 cell_tria->level_subdomain_id() != numbers::artificial_subdomain_id;
559
560 const typename DoFHandler<dim, spacedim>::level_cell_iterator cell_dofs(
561 &cell_tria->get_triangulation(),
562 cell_tria->level(),
563 cell_tria->index(),
564 &dof_handler);
565
566 const auto mapping_q =
567 dynamic_cast<const MappingQ<dim, spacedim> *>(&mapping);
568
569 // Step 2a) set up and reinit FEValues (if needed)
570 if (
571 ((vector_describes_relative_displacement ||
572 (is_non_artificial_cell == false)) &&
573 ((mapping_q != nullptr &&
574 this->get_degree() == mapping_q->get_degree()) ==
575 false)) /*condition 1: points need to be computed via FEValues*/
576 ||
577 (is_non_artificial_cell == true && interpolation_of_values_is_needed) /*condition 2: interpolation of values is needed*/)
578 {
579 // get FEValues (thread-safe); in the case that this thread has
580 // not created a an FEValues object yet, this helper-function also
581 // creates one with the right quadrature rule
582 auto &fe_values = fe_values_all.get();
583 if (fe_values.get() == nullptr)
584 {
585 QGaussLobatto<dim> quadrature_gl(this->polynomial_degree + 1);
586
587 std::vector<Point<dim>> quadrature_points;
588 for (const auto i :
589 FETools::hierarchic_to_lexicographic_numbering<dim>(
590 this->polynomial_degree))
591 quadrature_points.push_back(quadrature_gl.point(i));
592 Quadrature<dim> quadrature(quadrature_points);
593
594 fe_values = std::make_unique<FEValues<dim, spacedim>>(
595 mapping,
596 interpolation_of_values_is_needed ?
597 fe :
598 static_cast<const FiniteElement<dim, spacedim> &>(fe_nothing),
599 quadrature,
600 update_quadrature_points | update_values);
601 }
602
603 if (interpolation_of_values_is_needed)
604 fe_values->reinit(cell_dofs);
605 else
606 fe_values->reinit(cell_tria);
607 }
608
609 std::vector<Point<spacedim>> result;
610
611 // Step 2b) read of quadrature points in the relative displacement case
612 // note: we also take this path for non-active or artificial cells so that
613 // these cells are filled with some useful data
614 if (vector_describes_relative_displacement ||
615 (is_non_artificial_cell == false))
616 {
617 if (mapping_q != nullptr &&
618 this->get_degree() == mapping_q->get_degree())
619 result = mapping_q->compute_mapping_support_points(cell_tria);
620 else
621 result = fe_values_all.get()->get_quadrature_points();
622
623 // for non-active or artificial cells we are done here and return
624 // the absolute positions, since the provided vector cannot contain
625 // any useful information for these cells
626 if (is_non_artificial_cell == false)
627 return result;
628 }
629 else
630 {
631 result.resize(
632 Utilities::pow<unsigned int>(this->get_degree() + 1, dim));
633 }
634
635 // Step 2c) read global vector and adjust points accordingly
636 if (interpolation_of_values_is_needed == false)
637 {
638 // case 1: FE_Q or FE_DGQ with same degree as this class has; this
639 // is the simple case since no interpolation is needed
640 std::vector<types::global_dof_index> dof_indices(
641 fe.n_dofs_per_cell());
642 cell_dofs->get_mg_dof_indices(dof_indices);
643
644 for (unsigned int i = 0; i < dof_indices.size(); ++i)
645 {
646 const auto id = fe.system_to_component_index(i);
647
648 if (is_fe_q)
649 {
650 // case 1a: FE_Q
651 if (vector_describes_relative_displacement)
652 result[id.second][id.first] +=
653 vectors_ghosted[cell_tria->level()](dof_indices[i]);
654 else
655 result[id.second][id.first] =
656 vectors_ghosted[cell_tria->level()](dof_indices[i]);
657 }
658 else
659 {
660 // case 1b: FE_DGQ
661 if (vector_describes_relative_displacement)
662 result[lexicographic_to_hierarchic_numbering[id.second]]
663 [id.first] +=
664 vectors_ghosted[cell_tria->level()](dof_indices[i]);
665 else
666 result[lexicographic_to_hierarchic_numbering[id.second]]
667 [id.first] =
668 vectors_ghosted[cell_tria->level()](dof_indices[i]);
669 }
670 }
671 }
672 else
673 {
674 // case 2: general case; interpolation is needed
675 // note: the following code could be optimized for tensor-product
676 // elements via application of sum factorization as is done on
677 // MatrixFree/FEEvaluation
678 auto &fe_values = fe_values_all.get();
679
680 std::vector<types::global_dof_index> dof_indices(
681 fe.n_dofs_per_cell());
682 cell_dofs->get_mg_dof_indices(dof_indices);
683
684 std::vector<typename VectorType::value_type> dof_values(
685 fe.n_dofs_per_cell());
686
687 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
688 dof_values[i] = vectors_ghosted[cell_tria->level()](dof_indices[i]);
689
690 for (unsigned int c = 0; c < spacedim; ++c)
691 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
692 for (unsigned int q = 0; q < fe_values->n_quadrature_points; ++q)
693 if (vector_describes_relative_displacement == false && i == 0)
694 result[q][c] =
695 dof_values[i] * fe_values->shape_value_component(i, q, c);
696 else
697 result[q][c] +=
698 dof_values[i] * fe_values->shape_value_component(i, q, c);
699 }
700
701 return result;
702 });
703
704 uses_level_info = true;
705}
706
707
708
709template <int dim, int spacedim>
710std::size_t
712{
713 if (support_point_cache.get() != nullptr)
714 return sizeof(*this) +
715 MemoryConsumption::memory_consumption(*support_point_cache);
716 else
717 return sizeof(*this);
718}
719
720
721
722template <int dim, int spacedim>
723std::vector<Point<spacedim>>
725 const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
726{
727 Assert(support_point_cache.get() != nullptr,
728 ExcMessage("Must call MappingQCache::initialize() before "
729 "using it or after mesh has changed!"));
730
731 Assert(uses_level_info || cell->is_active(), ExcInternalError());
732
733 AssertIndexRange(cell->level(), support_point_cache->size());
734 AssertIndexRange(cell->index(), (*support_point_cache)[cell->level()].size());
735 return (*support_point_cache)[cell->level()][cell->index()];
736}
737
738
739
740template <int dim, int spacedim>
741boost::container::small_vector<Point<spacedim>,
744 const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
745{
746 Assert(support_point_cache.get() != nullptr,
747 ExcMessage("Must call MappingQCache::initialize() before "
748 "using it or after mesh has changed!"));
749
750 Assert(uses_level_info || cell->is_active(), ExcInternalError());
751
752 AssertIndexRange(cell->level(), support_point_cache->size());
753 AssertIndexRange(cell->index(), (*support_point_cache)[cell->level()].size());
754 const auto ptr = (*support_point_cache)[cell->level()][cell->index()].begin();
755 return boost::container::small_vector<Point<spacedim>,
757 ptr, ptr + cell->n_vertices());
758}
759
760
761
762//--------------------------- Explicit instantiations -----------------------
763#include "mapping_q_cache.inst"
764
765
const hp::FECollection< dim, spacedim > & get_fe_collection() const
const FiniteElement< dim, spacedim > & get_fe(const types::fe_index index=0) const
const IndexSet & locally_owned_mg_dofs(const unsigned int level) const
const Triangulation< dim, spacedim > & get_triangulation() const
const IndexSet & locally_owned_dofs() const
MPI_Comm get_communicator() const
Definition fe_q.h:551
const unsigned int degree
Definition fe_data.h:453
unsigned int n_dofs_per_cell() const
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
std::pair< unsigned int, unsigned int > system_to_component_index(const unsigned int index) const
unsigned int element_multiplicity(const unsigned int index) const
unsigned int n_base_elements() const
const unsigned int n_components
Definition function.h:164
virtual RangeNumberType value(const Point< dim > &p, const unsigned int component=0) const
void import_elements(const ::Vector< Number > &vec, VectorOperation::values operation, const std::shared_ptr< const Utilities::MPI::CommunicationPatternBase > &communication_pattern={})
virtual void reinit(const size_type size, const bool omit_zeroing_entries=false)
void import_elements(const Vector< Number, MemorySpace2 > &src, VectorOperation::values operation)
void reinit(const size_type size, const bool omit_zeroing_entries=false)
unsigned int max_level() const
unsigned int min_level() const
std::size_t memory_consumption() const
virtual boost::container::small_vector< Point< spacedim >, GeometryInfo< dim >::vertices_per_cell > get_vertices(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const override
virtual bool preserves_vertex_locations() const override
virtual std::vector< Point< spacedim > > compute_mapping_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const override
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
void initialize(const Mapping< dim, spacedim > &mapping, const Triangulation< dim, spacedim > &triangulation)
virtual std::vector< Point< spacedim > > compute_mapping_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const
friend class MappingQCache
Definition mapping_q.h:691
Abstract base class for mapping classes.
Definition mapping.h:317
static unsigned int n_threads()
Definition point.h:112
const Point< dim > & point(const unsigned int i) const
A class that provides a separate storage location on each thread that accesses the object.
virtual unsigned int n_global_levels() const
unsigned int size() const
Definition collection.h:265
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcMessage(std::string arg1)
@ update_quadrature_points
Transformed quadrature points.
IndexSet extract_locally_relevant_dofs(const DoFHandler< dim, spacedim > &dof_handler)
IndexSet extract_locally_relevant_level_dofs(const DoFHandler< dim, spacedim > &dof_handler, const unsigned int level)
std::enable_if_t< std::is_fundamental< T >::value, std::size_t > memory_consumption(const T &t)
std::vector< Integer > invert_permutation(const std::vector< Integer > &permutation)
Definition utilities.h:1672
void run(const std::vector< std::vector< Iterator > > &colored_iterators, Worker worker, Copier copier, const ScratchData &sample_scratch_data, const CopyData &sample_copy_data, const unsigned int queue_length=2 *MultithreadInfo::n_threads(), const unsigned int chunk_size=8)
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
const ::Triangulation< dim, spacedim > & tria