Loading [MathJax]/extensions/TeX/AMSsymbols.js
 deal.II version GIT relicensing-2669-g1eff4de72e 2025-02-20 01:30:00+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages Concepts
fe_values_views.cc
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 1998 - 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
17
19
20#include <deal.II/fe/fe.h>
24
25#include <deal.II/lac/vector.h>
26
27#ifdef DEAL_II_WITH_ADOLC
28# include <adolc/adouble.h>
29# include <adolc/adtl.h>
30#endif
31
32
34
35
36namespace internal
37{
38 namespace
39 {
40 template <int dim, int spacedim>
41 inline std::vector<unsigned int>
43 {
44 std::vector<unsigned int> shape_function_to_row_table(
45 fe.n_dofs_per_cell() * fe.n_components(),
47 unsigned int row = 0;
48 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
49 {
50 // loop over all components that are nonzero for this particular
51 // shape function. if a component is zero then we leave the
52 // value in the table unchanged (at the invalid value)
53 // otherwise it is mapped to the next free entry
54 unsigned int nth_nonzero_component = 0;
55 for (unsigned int c = 0; c < fe.n_components(); ++c)
56 if (fe.get_nonzero_components(i)[c] == true)
57 {
58 shape_function_to_row_table[i * fe.n_components() + c] =
59 row + nth_nonzero_component;
60 ++nth_nonzero_component;
61 }
62 row += fe.n_nonzero_components(i);
63 }
64
65 return shape_function_to_row_table;
66 }
67 } // namespace
68} // namespace internal
69
70
71
72namespace FEValuesViews
73{
74 template <int dim, int spacedim>
76 const unsigned int component)
77 : fe_values(&fe_values)
78 , component(component)
79 , shape_function_data(this->fe_values->fe->n_dofs_per_cell())
80 {
81 const FiniteElement<dim, spacedim> &fe = *this->fe_values->fe;
83
84 // TODO: we'd like to use the fields with the same name as these
85 // variables from FEValuesBase, but they aren't initialized yet
86 // at the time we get here, so re-create it all
87 const std::vector<unsigned int> shape_function_to_row_table =
89
90 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
91 {
92 const bool is_primitive = fe.is_primitive() || fe.is_primitive(i);
93
94 if (is_primitive == true)
95 shape_function_data[i].is_nonzero_shape_function_component =
96 (component == fe.system_to_component_index(i).first);
97 else
98 shape_function_data[i].is_nonzero_shape_function_component =
99 (fe.get_nonzero_components(i)[component] == true);
100
101 if (shape_function_data[i].is_nonzero_shape_function_component == true)
102 shape_function_data[i].row_index =
103 shape_function_to_row_table[i * fe.n_components() + component];
104 else
106 }
107 }
108
109
110
111 template <int dim, int spacedim>
113 : fe_values(nullptr)
114 , component(numbers::invalid_unsigned_int)
115 {}
116
117
118
119 template <int dim, int spacedim>
121 const unsigned int first_vector_component)
122 : fe_values(&fe_values)
123 , first_vector_component(first_vector_component)
124 , shape_function_data(this->fe_values->fe->n_dofs_per_cell())
125 {
126 const FiniteElement<dim, spacedim> &fe = *this->fe_values->fe;
128
129 // TODO: we'd like to use the fields with the same name as these
130 // variables from FEValuesBase, but they aren't initialized yet
131 // at the time we get here, so re-create it all
132 const std::vector<unsigned int> shape_function_to_row_table =
134
135 for (unsigned int d = 0; d < spacedim; ++d)
136 {
137 const unsigned int component = first_vector_component + d;
138
139 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
140 {
141 const bool is_primitive = fe.is_primitive() || fe.is_primitive(i);
142
143 if (is_primitive == true)
144 shape_function_data[i].is_nonzero_shape_function_component[d] =
145 (component == fe.system_to_component_index(i).first);
146 else
147 shape_function_data[i].is_nonzero_shape_function_component[d] =
148 (fe.get_nonzero_components(i)[component] == true);
149
150 if (shape_function_data[i].is_nonzero_shape_function_component[d] ==
151 true)
152 shape_function_data[i].row_index[d] =
153 shape_function_to_row_table[i * fe.n_components() + component];
154 else
155 shape_function_data[i].row_index[d] =
157 }
158 }
159
160 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
161 {
162 unsigned int n_nonzero_components = 0;
163 for (unsigned int d = 0; d < spacedim; ++d)
164 if (shape_function_data[i].is_nonzero_shape_function_component[d] ==
165 true)
166 ++n_nonzero_components;
167
168 if (n_nonzero_components == 0)
169 shape_function_data[i].single_nonzero_component = -2;
170 else if (n_nonzero_components > 1)
171 shape_function_data[i].single_nonzero_component = -1;
172 else
173 {
174 for (unsigned int d = 0; d < spacedim; ++d)
176 .is_nonzero_shape_function_component[d] == true)
177 {
178 shape_function_data[i].single_nonzero_component =
179 shape_function_data[i].row_index[d];
180 shape_function_data[i].single_nonzero_component_index = d;
181 break;
182 }
183 }
184 }
185 }
186
187
188
189 template <int dim, int spacedim>
191 : fe_values(nullptr)
192 , first_vector_component(numbers::invalid_unsigned_int)
193 {}
194
195
196
197 template <int dim, int spacedim>
199 const FEValuesBase<dim, spacedim> &fe_values,
200 const unsigned int first_tensor_component)
201 : fe_values(&fe_values)
202 , first_tensor_component(first_tensor_component)
203 , shape_function_data(this->fe_values->fe->n_dofs_per_cell())
204 {
205 const FiniteElement<dim, spacedim> &fe = *this->fe_values->fe;
206 Assert(first_tensor_component + (dim * dim + dim) / 2 - 1 <
207 fe.n_components(),
209 first_tensor_component +
211 0,
212 fe.n_components()));
213 // TODO: we'd like to use the fields with the same name as these
214 // variables from FEValuesBase, but they aren't initialized yet
215 // at the time we get here, so re-create it all
216 const std::vector<unsigned int> shape_function_to_row_table =
218
219 for (unsigned int d = 0;
220 d < ::SymmetricTensor<2, dim>::n_independent_components;
221 ++d)
222 {
223 const unsigned int component = first_tensor_component + d;
224
225 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
226 {
227 const bool is_primitive = fe.is_primitive() || fe.is_primitive(i);
228
229 if (is_primitive == true)
230 shape_function_data[i].is_nonzero_shape_function_component[d] =
231 (component == fe.system_to_component_index(i).first);
232 else
233 shape_function_data[i].is_nonzero_shape_function_component[d] =
234 (fe.get_nonzero_components(i)[component] == true);
235
236 if (shape_function_data[i].is_nonzero_shape_function_component[d] ==
237 true)
238 shape_function_data[i].row_index[d] =
239 shape_function_to_row_table[i * fe.n_components() + component];
240 else
241 shape_function_data[i].row_index[d] =
243 }
244 }
245
246 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
247 {
248 unsigned int n_nonzero_components = 0;
249 for (unsigned int d = 0;
250 d < ::SymmetricTensor<2, dim>::n_independent_components;
251 ++d)
252 if (shape_function_data[i].is_nonzero_shape_function_component[d] ==
253 true)
254 ++n_nonzero_components;
255
256 if (n_nonzero_components == 0)
257 shape_function_data[i].single_nonzero_component = -2;
258 else if (n_nonzero_components > 1)
259 shape_function_data[i].single_nonzero_component = -1;
260 else
261 {
262 for (unsigned int d = 0;
263 d < ::SymmetricTensor<2, dim>::n_independent_components;
264 ++d)
265 if (shape_function_data[i]
266 .is_nonzero_shape_function_component[d] == true)
267 {
268 shape_function_data[i].single_nonzero_component =
269 shape_function_data[i].row_index[d];
270 shape_function_data[i].single_nonzero_component_index = d;
271 break;
272 }
273 }
274 }
275 }
276
277
278
279 template <int dim, int spacedim>
281 : fe_values(nullptr)
282 , first_tensor_component(numbers::invalid_unsigned_int)
283 {}
284
285
286
287 template <int dim, int spacedim>
289 const unsigned int first_tensor_component)
290 : fe_values(&fe_values)
291 , first_tensor_component(first_tensor_component)
292 , shape_function_data(this->fe_values->fe->n_dofs_per_cell())
293 {
294 const FiniteElement<dim, spacedim> &fe = *this->fe_values->fe;
295 AssertIndexRange(first_tensor_component + dim * dim - 1, fe.n_components());
296 // TODO: we'd like to use the fields with the same name as these
297 // variables from FEValuesBase, but they aren't initialized yet
298 // at the time we get here, so re-create it all
299 const std::vector<unsigned int> shape_function_to_row_table =
301
302 for (unsigned int d = 0; d < dim * dim; ++d)
303 {
304 const unsigned int component = first_tensor_component + d;
305
306 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
307 {
308 const bool is_primitive = fe.is_primitive() || fe.is_primitive(i);
309
310 if (is_primitive == true)
311 shape_function_data[i].is_nonzero_shape_function_component[d] =
312 (component == fe.system_to_component_index(i).first);
313 else
314 shape_function_data[i].is_nonzero_shape_function_component[d] =
315 (fe.get_nonzero_components(i)[component] == true);
316
317 if (shape_function_data[i].is_nonzero_shape_function_component[d] ==
318 true)
319 shape_function_data[i].row_index[d] =
320 shape_function_to_row_table[i * fe.n_components() + component];
321 else
322 shape_function_data[i].row_index[d] =
324 }
325 }
326
327 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
328 {
329 unsigned int n_nonzero_components = 0;
330 for (unsigned int d = 0; d < dim * dim; ++d)
331 if (shape_function_data[i].is_nonzero_shape_function_component[d] ==
332 true)
333 ++n_nonzero_components;
334
335 if (n_nonzero_components == 0)
336 shape_function_data[i].single_nonzero_component = -2;
337 else if (n_nonzero_components > 1)
338 shape_function_data[i].single_nonzero_component = -1;
339 else
340 {
341 for (unsigned int d = 0; d < dim * dim; ++d)
342 if (shape_function_data[i]
343 .is_nonzero_shape_function_component[d] == true)
344 {
345 shape_function_data[i].single_nonzero_component =
346 shape_function_data[i].row_index[d];
347 shape_function_data[i].single_nonzero_component_index = d;
348 break;
349 }
350 }
351 }
352 }
353
354
355
356 template <int dim, int spacedim>
358 : fe_values(nullptr)
359 , first_tensor_component(numbers::invalid_unsigned_int)
360 {}
361
362
363
364 template <int dim, int spacedim>
365 template <typename Number>
366 void
368 const ReadVector<Number> &fe_function,
369 std::vector<solution_value_type<Number>> &values) const
370 {
371 Assert(fe_values->update_flags & update_values,
373 "update_values")));
374 Assert(fe_values->present_cell.is_initialized(),
376 AssertDimension(values.size(), fe_values->n_quadrature_points);
377
378 // get function values of dofs on this cell and call internal worker
379 // function
380 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
381 fe_values->present_cell.get_interpolated_dof_values(fe_function,
382 dof_values);
383 internal::do_function_values<dim, spacedim>(
384 make_const_array_view(dof_values),
385 fe_values->finite_element_output.shape_values,
386 shape_function_data,
387 values);
388 }
389
390
391
392 template <int dim, int spacedim>
393 template <class InputVector>
394 void
396 const InputVector &dof_values,
398 const
399 {
400 Assert(fe_values->update_flags & update_values,
402 "update_values")));
403 Assert(fe_values->present_cell.is_initialized(),
405 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
406 AssertDimension(values.size(), fe_values->n_quadrature_points);
407
408 internal::do_function_values<dim, spacedim>(
409 make_const_array_view(dof_values),
410 fe_values->finite_element_output.shape_values,
411 shape_function_data,
412 values);
413 }
414
415
416
417 template <int dim, int spacedim>
418 template <typename Number>
419 void
421 const ReadVector<Number> &fe_function,
422 std::vector<solution_gradient_type<Number>> &gradients) const
423 {
424 Assert(fe_values->update_flags & update_gradients,
426 "update_gradients")));
427 Assert(fe_values->present_cell.is_initialized(),
429 AssertDimension(fe_function.size(),
430 fe_values->present_cell.n_dofs_for_dof_handler());
431 AssertDimension(gradients.size(), fe_values->n_quadrature_points);
432
433 // get function values of dofs on this cell
434 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
435 fe_values->present_cell.get_interpolated_dof_values(fe_function,
436 dof_values);
437 internal::do_function_derivatives<1, dim, spacedim>(
438 make_const_array_view(dof_values),
439 fe_values->finite_element_output.shape_gradients,
440 shape_function_data,
441 gradients);
442 }
443
444
445
446 template <int dim, int spacedim>
447 template <typename InputVector>
448 void
450 const InputVector &dof_values,
452 &gradients) const
453 {
454 Assert(fe_values->update_flags & update_gradients,
456 "update_gradients")));
457 Assert(fe_values->present_cell.is_initialized(),
459 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
460 AssertDimension(gradients.size(), fe_values->n_quadrature_points);
461
462 internal::do_function_derivatives<1, dim, spacedim>(
463 make_const_array_view(dof_values),
464 fe_values->finite_element_output.shape_gradients,
465 shape_function_data,
466 gradients);
467 }
468
469
470
471 template <int dim, int spacedim>
472 template <typename Number>
473 void
475 const ReadVector<Number> &fe_function,
476 std::vector<solution_hessian_type<Number>> &hessians) const
477 {
478 Assert(fe_values->update_flags & update_hessians,
480 "update_hessians")));
481 Assert(fe_values->present_cell.is_initialized(),
483 AssertDimension(fe_function.size(),
484 fe_values->present_cell.n_dofs_for_dof_handler());
485 AssertDimension(hessians.size(), fe_values->n_quadrature_points);
486
487 // get function values of dofs on this cell
488 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
489 fe_values->present_cell.get_interpolated_dof_values(fe_function,
490 dof_values);
491 internal::do_function_derivatives<2, dim, spacedim>(
492 make_const_array_view(dof_values),
493 fe_values->finite_element_output.shape_hessians,
494 shape_function_data,
495 hessians);
496 }
497
498
499
500 template <int dim, int spacedim>
501 template <class InputVector>
502 void
504 const InputVector &dof_values,
506 &hessians) const
507 {
508 Assert(fe_values->update_flags & update_hessians,
510 "update_hessians")));
511 Assert(fe_values->present_cell.is_initialized(),
513 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
514 AssertDimension(hessians.size(), fe_values->n_quadrature_points);
515
516 internal::do_function_derivatives<2, dim, spacedim>(
517 make_const_array_view(dof_values),
518 fe_values->finite_element_output.shape_hessians,
519 shape_function_data,
520 hessians);
521 }
522
523
524
525 template <int dim, int spacedim>
526 template <typename Number>
527 void
529 const ReadVector<Number> &fe_function,
530 std::vector<solution_laplacian_type<Number>> &laplacians) const
531 {
532 Assert(fe_values->update_flags & update_hessians,
534 "update_hessians")));
535 Assert(fe_values->present_cell.is_initialized(),
537 AssertDimension(fe_function.size(),
538 fe_values->present_cell.n_dofs_for_dof_handler());
539 AssertDimension(laplacians.size(), fe_values->n_quadrature_points);
540
541 // get function values of dofs on this cell
542 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
543 fe_values->present_cell.get_interpolated_dof_values(fe_function,
544 dof_values);
545 internal::do_function_laplacians<dim, spacedim>(
546 make_const_array_view(dof_values),
547 fe_values->finite_element_output.shape_hessians,
548 shape_function_data,
549 laplacians);
550 }
551
552
553
554 template <int dim, int spacedim>
555 template <class InputVector>
556 void
558 const InputVector &dof_values,
560 &laplacians) const
561 {
562 Assert(fe_values->update_flags & update_hessians,
564 "update_hessians")));
565 Assert(fe_values->present_cell.is_initialized(),
567 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
568 AssertDimension(laplacians.size(), fe_values->n_quadrature_points);
569
570 internal::do_function_laplacians<dim, spacedim>(
571 make_const_array_view(dof_values),
572 fe_values->finite_element_output.shape_hessians,
573 shape_function_data,
574 laplacians);
575 }
576
577
578
579 template <int dim, int spacedim>
580 template <typename Number>
581 void
583 const ReadVector<Number> &fe_function,
584 std::vector<solution_third_derivative_type<Number>> &third_derivatives)
585 const
586 {
587 Assert(fe_values->update_flags & update_3rd_derivatives,
589 "update_3rd_derivatives")));
590 Assert(fe_values->present_cell.is_initialized(),
592 AssertDimension(fe_function.size(),
593 fe_values->present_cell.n_dofs_for_dof_handler());
594 AssertDimension(third_derivatives.size(), fe_values->n_quadrature_points);
595
596 // get function values of dofs on this cell
597 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
598 fe_values->present_cell.get_interpolated_dof_values(fe_function,
599 dof_values);
600 internal::do_function_derivatives<3, dim, spacedim>(
601 make_const_array_view(dof_values),
602 fe_values->finite_element_output.shape_3rd_derivatives,
603 shape_function_data,
604 third_derivatives);
605 }
606
607
608
609 template <int dim, int spacedim>
610 template <class InputVector>
611 void
613 const InputVector &dof_values,
614 std::vector<
616 &third_derivatives) const
617 {
618 Assert(fe_values->update_flags & update_3rd_derivatives,
620 "update_3rd_derivatives")));
621 Assert(fe_values->present_cell.is_initialized(),
623 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
624 AssertDimension(third_derivatives.size(), fe_values->n_quadrature_points);
625
626 internal::do_function_derivatives<3, dim, spacedim>(
627 make_const_array_view(dof_values),
628 fe_values->finite_element_output.shape_3rd_derivatives,
629 shape_function_data,
630 third_derivatives);
631 }
632
633
634
635 template <int dim, int spacedim>
636 template <typename Number>
637 void
639 const ReadVector<Number> &fe_function,
640 std::vector<solution_value_type<Number>> &values) const
641 {
642 Assert(fe_values->update_flags & update_values,
644 "update_values")));
645 Assert(fe_values->present_cell.is_initialized(),
647 AssertDimension(values.size(), fe_values->n_quadrature_points);
648
649 // get function values of dofs on this cell
650 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
651 fe_values->present_cell.get_interpolated_dof_values(fe_function,
652 dof_values);
653 internal::do_function_values<dim, spacedim>(
654 make_const_array_view(dof_values),
655 fe_values->finite_element_output.shape_values,
656 shape_function_data,
657 values);
658 }
659
660
661
662 template <int dim, int spacedim>
663 template <class InputVector>
664 void
666 const InputVector &dof_values,
668 const
669 {
670 Assert(fe_values->update_flags & update_values,
672 "update_values")));
673 Assert(fe_values->present_cell.is_initialized(),
675 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
676 AssertDimension(values.size(), fe_values->n_quadrature_points);
677
678 internal::do_function_values<dim, spacedim>(
679 make_const_array_view(dof_values),
680 fe_values->finite_element_output.shape_values,
681 shape_function_data,
682 values);
683 }
684
685
686
687 template <int dim, int spacedim>
688 template <typename Number>
689 void
691 const ReadVector<Number> &fe_function,
692 std::vector<solution_gradient_type<Number>> &gradients) const
693 {
694 Assert(fe_values->update_flags & update_gradients,
696 "update_gradients")));
697 Assert(fe_values->present_cell.is_initialized(),
699 AssertDimension(fe_function.size(),
700 fe_values->present_cell.n_dofs_for_dof_handler());
701 AssertDimension(gradients.size(), fe_values->n_quadrature_points);
702
703 // get function values of dofs on this cell
704 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
705 fe_values->present_cell.get_interpolated_dof_values(fe_function,
706 dof_values);
707 internal::do_function_derivatives<1, dim, spacedim>(
708 make_const_array_view(dof_values),
709 fe_values->finite_element_output.shape_gradients,
710 shape_function_data,
711 gradients);
712 }
713
714
715
716 template <int dim, int spacedim>
717 template <typename InputVector>
718 void
720 const InputVector &dof_values,
722 &gradients) const
723 {
724 Assert(fe_values->update_flags & update_gradients,
726 "update_gradients")));
727 Assert(fe_values->present_cell.is_initialized(),
729 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
730 AssertDimension(gradients.size(), fe_values->n_quadrature_points);
731
732 internal::do_function_derivatives<1, dim, spacedim>(
733 make_const_array_view(dof_values),
734 fe_values->finite_element_output.shape_gradients,
735 shape_function_data,
736 gradients);
737 }
738
739
740
741 template <int dim, int spacedim>
742 template <typename Number>
743 void
745 const ReadVector<Number> &fe_function,
746 std::vector<solution_symmetric_gradient_type<Number>> &symmetric_gradients)
747 const
748 {
749 Assert(fe_values->update_flags & update_gradients,
751 "update_gradients")));
752 Assert(fe_values->present_cell.is_initialized(),
754 AssertDimension(fe_function.size(),
755 fe_values->present_cell.n_dofs_for_dof_handler());
756 AssertDimension(symmetric_gradients.size(), fe_values->n_quadrature_points);
757
758 // get function values of dofs on this cell
759 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
760 fe_values->present_cell.get_interpolated_dof_values(fe_function,
761 dof_values);
762 internal::do_function_symmetric_gradients<dim, spacedim>(
763 make_const_array_view(dof_values),
764 fe_values->finite_element_output.shape_gradients,
765 shape_function_data,
766 symmetric_gradients);
767 }
768
769
770
771 template <int dim, int spacedim>
772 template <class InputVector>
773 void
775 const InputVector &dof_values,
776 std::vector<
778 &symmetric_gradients) const
779 {
780 Assert(fe_values->update_flags & update_gradients,
782 "update_gradients")));
783 Assert(fe_values->present_cell.is_initialized(),
785 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
786 AssertDimension(symmetric_gradients.size(), fe_values->n_quadrature_points);
787
788 internal::do_function_symmetric_gradients<dim, spacedim>(
789 make_const_array_view(dof_values),
790 fe_values->finite_element_output.shape_gradients,
791 shape_function_data,
792 symmetric_gradients);
793 }
794
795
796
797 template <int dim, int spacedim>
798 template <typename Number>
799 void
801 const ReadVector<Number> &fe_function,
802 std::vector<solution_divergence_type<Number>> &divergences) const
803 {
804 Assert(fe_values->update_flags & update_gradients,
806 "update_gradients")));
807 Assert(fe_values->present_cell.is_initialized(),
809 AssertDimension(fe_function.size(),
810 fe_values->present_cell.n_dofs_for_dof_handler());
811 AssertDimension(divergences.size(), fe_values->n_quadrature_points);
812
813 // get function values of dofs
814 // on this cell
815 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
816 fe_values->present_cell.get_interpolated_dof_values(fe_function,
817 dof_values);
818 internal::do_function_divergences<dim, spacedim>(
819 make_const_array_view(dof_values),
820 fe_values->finite_element_output.shape_gradients,
821 shape_function_data,
822 divergences);
823 }
824
825
826
827 template <int dim, int spacedim>
828 template <class InputVector>
829 void
831 const InputVector &dof_values,
833 &divergences) const
834 {
835 Assert(fe_values->update_flags & update_gradients,
837 "update_gradients")));
838 Assert(fe_values->present_cell.is_initialized(),
840 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
841 AssertDimension(divergences.size(), fe_values->n_quadrature_points);
842
843 internal::do_function_divergences<dim, spacedim>(
844 make_const_array_view(dof_values),
845 fe_values->finite_element_output.shape_gradients,
846 shape_function_data,
847 divergences);
848 }
849
850
851
852 template <int dim, int spacedim>
853 template <typename Number>
854 void
856 const ReadVector<Number> &fe_function,
857 std::vector<solution_curl_type<Number>> &curls) const
858 {
859 Assert(fe_values->update_flags & update_gradients,
861 "update_gradients")));
862 Assert(fe_values->present_cell.is_initialized(),
863 ExcMessage("FEValues object is not reinited to any cell"));
864 AssertDimension(fe_function.size(),
865 fe_values->present_cell.n_dofs_for_dof_handler());
866 AssertDimension(curls.size(), fe_values->n_quadrature_points);
867
868 // get function values of dofs on this cell
869 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
870 fe_values->present_cell.get_interpolated_dof_values(fe_function,
871 dof_values);
872 internal::do_function_curls<dim, spacedim>(
873 make_const_array_view(dof_values),
874 fe_values->finite_element_output.shape_gradients,
875 shape_function_data,
876 curls);
877 }
878
879
880
881 template <int dim, int spacedim>
882 template <class InputVector>
883 void
885 const InputVector &dof_values,
887 const
888 {
889 Assert(fe_values->update_flags & update_gradients,
891 "update_gradients")));
892 Assert(fe_values->present_cell.is_initialized(),
893 ExcMessage("FEValues object is not reinited to any cell"));
894 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
895 AssertDimension(curls.size(), fe_values->n_quadrature_points);
896
897 internal::do_function_curls<dim, spacedim>(
898 make_const_array_view(dof_values),
899 fe_values->finite_element_output.shape_gradients,
900 shape_function_data,
901 curls);
902 }
903
904
905
906 template <int dim, int spacedim>
907 template <typename Number>
908 void
910 const ReadVector<Number> &fe_function,
911 std::vector<solution_hessian_type<Number>> &hessians) const
912 {
913 Assert(fe_values->update_flags & update_hessians,
915 "update_hessians")));
916 Assert(fe_values->present_cell.is_initialized(),
918 AssertDimension(fe_function.size(),
919 fe_values->present_cell.n_dofs_for_dof_handler());
920 AssertDimension(hessians.size(), fe_values->n_quadrature_points);
921
922 // get function values of dofs on this cell
923 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
924 fe_values->present_cell.get_interpolated_dof_values(fe_function,
925 dof_values);
926 internal::do_function_derivatives<2, dim, spacedim>(
927 make_const_array_view(dof_values),
928 fe_values->finite_element_output.shape_hessians,
929 shape_function_data,
930 hessians);
931 }
932
933
934
935 template <int dim, int spacedim>
936 template <class InputVector>
937 void
939 const InputVector &dof_values,
941 &hessians) const
942 {
943 Assert(fe_values->update_flags & update_hessians,
945 "update_hessians")));
946 Assert(fe_values->present_cell.is_initialized(),
948 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
949 AssertDimension(hessians.size(), fe_values->n_quadrature_points);
950
951 internal::do_function_derivatives<2, dim, spacedim>(
952 make_const_array_view(dof_values),
953 fe_values->finite_element_output.shape_hessians,
954 shape_function_data,
955 hessians);
956 }
957
958
959
960 template <int dim, int spacedim>
961 template <typename Number>
962 void
964 const ReadVector<Number> &fe_function,
965 std::vector<solution_value_type<Number>> &laplacians) const
966 {
967 Assert(fe_values->update_flags & update_hessians,
969 "update_hessians")));
970 Assert(laplacians.size() == fe_values->n_quadrature_points,
971 ExcDimensionMismatch(laplacians.size(),
972 fe_values->n_quadrature_points));
973 Assert(fe_values->present_cell.is_initialized(),
975 Assert(
976 fe_function.size() == fe_values->present_cell.n_dofs_for_dof_handler(),
977 ExcDimensionMismatch(fe_function.size(),
978 fe_values->present_cell.n_dofs_for_dof_handler()));
979 AssertDimension(laplacians.size(), fe_values->n_quadrature_points);
980
981 // get function values of dofs on this cell
982 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
983 fe_values->present_cell.get_interpolated_dof_values(fe_function,
984 dof_values);
985 internal::do_function_laplacians<dim, spacedim>(
986 make_const_array_view(dof_values),
987 fe_values->finite_element_output.shape_hessians,
988 shape_function_data,
989 laplacians);
990 }
991
992
993
994 template <int dim, int spacedim>
995 template <class InputVector>
996 void
998 const InputVector &dof_values,
1000 &laplacians) const
1001 {
1002 Assert(fe_values->update_flags & update_hessians,
1004 "update_hessians")));
1005 Assert(laplacians.size() == fe_values->n_quadrature_points,
1006 ExcDimensionMismatch(laplacians.size(),
1007 fe_values->n_quadrature_points));
1008 Assert(fe_values->present_cell.is_initialized(),
1010 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
1011 AssertDimension(laplacians.size(), fe_values->n_quadrature_points);
1012
1013 internal::do_function_laplacians<dim, spacedim>(
1014 make_const_array_view(dof_values),
1015 fe_values->finite_element_output.shape_hessians,
1016 shape_function_data,
1017 laplacians);
1018 }
1019
1020
1021
1022 template <int dim, int spacedim>
1023 template <typename Number>
1024 void
1026 const ReadVector<Number> &fe_function,
1027 std::vector<solution_third_derivative_type<Number>> &third_derivatives)
1028 const
1029 {
1030 Assert(fe_values->update_flags & update_3rd_derivatives,
1032 "update_3rd_derivatives")));
1033 Assert(fe_values->present_cell.is_initialized(),
1035 AssertDimension(fe_function.size(),
1036 fe_values->present_cell.n_dofs_for_dof_handler());
1037 AssertDimension(third_derivatives.size(), fe_values->n_quadrature_points);
1038
1039 // get function values of dofs on this cell
1040 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
1041 fe_values->present_cell.get_interpolated_dof_values(fe_function,
1042 dof_values);
1043 internal::do_function_derivatives<3, dim, spacedim>(
1044 make_const_array_view(dof_values),
1045 fe_values->finite_element_output.shape_3rd_derivatives,
1046 shape_function_data,
1047 third_derivatives);
1048 }
1049
1050
1051
1052 template <int dim, int spacedim>
1053 template <class InputVector>
1054 void
1056 const InputVector &dof_values,
1057 std::vector<
1059 &third_derivatives) const
1060 {
1061 Assert(fe_values->update_flags & update_3rd_derivatives,
1063 "update_3rd_derivatives")));
1064 Assert(fe_values->present_cell.is_initialized(),
1066 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
1067 AssertDimension(third_derivatives.size(), fe_values->n_quadrature_points);
1068
1069 internal::do_function_derivatives<3, dim, spacedim>(
1070 make_const_array_view(dof_values),
1071 fe_values->finite_element_output.shape_3rd_derivatives,
1072 shape_function_data,
1073 third_derivatives);
1074 }
1075
1076
1077
1078 template <int dim, int spacedim>
1079 template <typename Number>
1080 void
1082 const ReadVector<Number> &fe_function,
1083 std::vector<solution_value_type<Number>> &values) const
1084 {
1085 Assert(fe_values->update_flags & update_values,
1087 "update_values")));
1088 Assert(fe_values->present_cell.is_initialized(),
1090 AssertDimension(values.size(), fe_values->n_quadrature_points);
1091
1092 // get function values of dofs on this cell
1093 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
1094 fe_values->present_cell.get_interpolated_dof_values(fe_function,
1095 dof_values);
1096 internal::do_function_values<dim, spacedim>(
1097 make_const_array_view(dof_values),
1098 fe_values->finite_element_output.shape_values,
1099 shape_function_data,
1100 values);
1101 }
1102
1103
1104
1105 template <int dim, int spacedim>
1106 template <class InputVector>
1107 void
1109 const InputVector &dof_values,
1111 const
1112 {
1113 Assert(fe_values->update_flags & update_values,
1115 "update_values")));
1116 Assert(fe_values->present_cell.is_initialized(),
1118 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
1119 AssertDimension(values.size(), fe_values->n_quadrature_points);
1120
1121 internal::do_function_values<dim, spacedim>(
1122 make_const_array_view(dof_values),
1123 fe_values->finite_element_output.shape_values,
1124 shape_function_data,
1125 values);
1126 }
1127
1128
1129
1130 template <int dim, int spacedim>
1131 template <typename Number>
1132 void
1134 const ReadVector<Number> &fe_function,
1135 std::vector<solution_divergence_type<Number>> &divergences) const
1136 {
1137 Assert(fe_values->update_flags & update_gradients,
1139 "update_gradients")));
1140 Assert(fe_values->present_cell.is_initialized(),
1142 AssertDimension(fe_function.size(),
1143 fe_values->present_cell.n_dofs_for_dof_handler());
1144 AssertDimension(divergences.size(), fe_values->n_quadrature_points);
1145
1146 // get function values of dofs
1147 // on this cell
1148 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
1149 fe_values->present_cell.get_interpolated_dof_values(fe_function,
1150 dof_values);
1151 internal::do_function_divergences<dim, spacedim>(
1152 make_const_array_view(dof_values),
1153 fe_values->finite_element_output.shape_gradients,
1154 shape_function_data,
1155 divergences);
1156 }
1157
1158
1159
1160 template <int dim, int spacedim>
1161 template <class InputVector>
1162 void
1165 const InputVector &dof_values,
1167 &divergences) const
1168 {
1169 Assert(fe_values->update_flags & update_gradients,
1171 "update_gradients")));
1172 Assert(fe_values->present_cell.is_initialized(),
1174 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
1175 AssertDimension(divergences.size(), fe_values->n_quadrature_points);
1176
1177 internal::do_function_divergences<dim, spacedim>(
1178 make_const_array_view(dof_values),
1179 fe_values->finite_element_output.shape_gradients,
1180 shape_function_data,
1181 divergences);
1182 }
1183
1184
1185
1186 template <int dim, int spacedim>
1187 template <typename Number>
1188 void
1190 const ReadVector<Number> &fe_function,
1191 std::vector<solution_value_type<Number>> &values) const
1192 {
1193 Assert(fe_values->update_flags & update_values,
1195 "update_values")));
1196 Assert(fe_values->present_cell.is_initialized(),
1198 AssertDimension(values.size(), fe_values->n_quadrature_points);
1199
1200 // get function values of dofs on this cell
1201 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
1202 fe_values->present_cell.get_interpolated_dof_values(fe_function,
1203 dof_values);
1204 internal::do_function_values<dim, spacedim>(
1205 make_const_array_view(dof_values),
1206 fe_values->finite_element_output.shape_values,
1207 shape_function_data,
1208 values);
1209 }
1210
1211
1212
1213 template <int dim, int spacedim>
1214 template <class InputVector>
1215 void
1217 const InputVector &dof_values,
1219 const
1220 {
1221 Assert(fe_values->update_flags & update_values,
1223 "update_values")));
1224 Assert(fe_values->present_cell.is_initialized(),
1226 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
1227 AssertDimension(values.size(), fe_values->n_quadrature_points);
1228
1229 internal::do_function_values<dim, spacedim>(
1230 make_const_array_view(dof_values),
1231 fe_values->finite_element_output.shape_values,
1232 shape_function_data,
1233 values);
1234 }
1235
1236
1237
1238 template <int dim, int spacedim>
1239 template <typename Number>
1240 void
1242 const ReadVector<Number> &fe_function,
1243 std::vector<solution_divergence_type<Number>> &divergences) const
1244 {
1245 Assert(fe_values->update_flags & update_gradients,
1247 "update_gradients")));
1248 Assert(fe_values->present_cell.is_initialized(),
1250 AssertDimension(fe_function.size(),
1251 fe_values->present_cell.n_dofs_for_dof_handler());
1252 AssertDimension(divergences.size(), fe_values->n_quadrature_points);
1253
1254 // get function values of dofs
1255 // on this cell
1256 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
1257 fe_values->present_cell.get_interpolated_dof_values(fe_function,
1258 dof_values);
1259 internal::do_function_divergences<dim, spacedim>(
1260 make_const_array_view(dof_values),
1261 fe_values->finite_element_output.shape_gradients,
1262 shape_function_data,
1263 divergences);
1264 }
1265
1266
1267
1268 template <int dim, int spacedim>
1269 template <class InputVector>
1270 void
1272 const InputVector &dof_values,
1274 &divergences) const
1275 {
1276 Assert(fe_values->update_flags & update_gradients,
1278 "update_gradients")));
1279 Assert(fe_values->present_cell.is_initialized(),
1281 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
1282 AssertDimension(divergences.size(), fe_values->n_quadrature_points);
1283
1284 internal::do_function_divergences<dim, spacedim>(
1285 make_const_array_view(dof_values),
1286 fe_values->finite_element_output.shape_gradients,
1287 shape_function_data,
1288 divergences);
1289 }
1290
1291
1292
1293 template <int dim, int spacedim>
1294 template <typename Number>
1295 void
1297 const ReadVector<Number> &fe_function,
1298 std::vector<solution_gradient_type<Number>> &gradients) const
1299 {
1300 Assert(fe_values->update_flags & update_gradients,
1302 "update_gradients")));
1303 Assert(fe_values->present_cell.is_initialized(),
1305 AssertDimension(fe_function.size(),
1306 fe_values->present_cell.n_dofs_for_dof_handler());
1307 AssertDimension(gradients.size(), fe_values->n_quadrature_points);
1308
1309 // get function values of dofs
1310 // on this cell
1311 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
1312 fe_values->present_cell.get_interpolated_dof_values(fe_function,
1313 dof_values);
1314 internal::do_function_gradients<dim, spacedim>(
1315 make_const_array_view(dof_values),
1316 fe_values->finite_element_output.shape_gradients,
1317 shape_function_data,
1318 gradients);
1319 }
1320
1321
1322
1323 template <int dim, int spacedim>
1324 template <class InputVector>
1325 void
1327 const InputVector &dof_values,
1329 &gradients) const
1330 {
1331 Assert(fe_values->update_flags & update_gradients,
1333 "update_gradients")));
1334 Assert(fe_values->present_cell.is_initialized(),
1336 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
1337 AssertDimension(gradients.size(), fe_values->n_quadrature_points);
1338
1339 internal::do_function_gradients<dim, spacedim>(
1340 make_const_array_view(dof_values),
1341 fe_values->finite_element_output.shape_gradients,
1342 shape_function_data,
1343 gradients);
1344 }
1345} // namespace FEValuesViews
1346
1347
1348namespace internal
1349{
1350 namespace FEValuesViews
1351 {
1352 template <int dim, int spacedim>
1354 {
1355 const FiniteElement<dim, spacedim> &fe = fe_values.get_fe();
1356
1357 const unsigned int n_scalars = fe.n_components();
1358 scalars.resize(n_scalars);
1359
1360 // compute number of vectors that we can fit into this finite element.
1361 // note that this is based on the dimensionality 'dim' of the manifold,
1362 // not 'spacedim' of the output vector
1363 const unsigned int n_vectors =
1366 1 :
1367 0);
1368 vectors.resize(n_vectors);
1369
1370 // compute number of symmetric tensors in the same way as above
1371 const unsigned int n_symmetric_second_order_tensors =
1372 (fe.n_components() >=
1374 fe.n_components() -
1376 0);
1377 symmetric_second_order_tensors.resize(n_symmetric_second_order_tensors);
1378
1379 // compute number of symmetric tensors in the same way as above
1380 const unsigned int n_second_order_tensors =
1383 1 :
1384 0);
1385 second_order_tensors.resize(n_second_order_tensors);
1386 }
1387 } // namespace FEValuesViews
1388} // namespace internal
1389
1390/*------------------------------- Explicit Instantiations -------------*/
1391
1392#include "fe/fe_values_views.inst"
1393
auto make_const_array_view(const Container &container) -> decltype(make_array_view(container))
const ObserverPointer< const FiniteElement< dim, spacedim >, FEValuesBase< dim, spacedim > > fe
const FiniteElement< dim, spacedim > & get_fe() const
void get_function_third_derivatives_from_local_dof_values(const InputVector &dof_values, std::vector< solution_third_derivative_type< typename InputVector::value_type > > &third_derivatives) const
void get_function_laplacians(const ReadVector< Number > &fe_function, std::vector< solution_laplacian_type< Number > > &laplacians) const
typename ProductType< Number, hessian_type >::type solution_hessian_type
void get_function_hessians_from_local_dof_values(const InputVector &dof_values, std::vector< solution_hessian_type< typename InputVector::value_type > > &hessians) const
void get_function_values_from_local_dof_values(const InputVector &dof_values, std::vector< solution_value_type< typename InputVector::value_type > > &values) const
std::vector< ShapeFunctionData > shape_function_data
void get_function_gradients(const ReadVector< Number > &fe_function, std::vector< solution_gradient_type< Number > > &gradients) const
typename ProductType< Number, value_type >::type solution_laplacian_type
typename ProductType< Number, third_derivative_type >::type solution_third_derivative_type
typename ProductType< Number, value_type >::type solution_value_type
typename ProductType< Number, gradient_type >::type solution_gradient_type
void get_function_values(const ReadVector< Number > &fe_function, std::vector< solution_value_type< Number > > &values) const
void get_function_laplacians_from_local_dof_values(const InputVector &dof_values, std::vector< solution_laplacian_type< typename InputVector::value_type > > &laplacians) const
void get_function_gradients_from_local_dof_values(const InputVector &dof_values, std::vector< solution_gradient_type< typename InputVector::value_type > > &gradients) const
void get_function_third_derivatives(const ReadVector< Number > &fe_function, std::vector< solution_third_derivative_type< Number > > &third_derivatives) const
void get_function_hessians(const ReadVector< Number > &fe_function, std::vector< solution_hessian_type< Number > > &hessians) const
typename ProductType< Number, value_type >::type solution_value_type
typename ProductType< Number, divergence_type >::type solution_divergence_type
typename ProductType< Number, value_type >::type solution_value_type
typename ProductType< Number, divergence_type >::type solution_divergence_type
typename ProductType< Number, gradient_type >::type solution_gradient_type
typename ProductType< Number, third_derivative_type >::type solution_third_derivative_type
void get_function_symmetric_gradients_from_local_dof_values(const InputVector &dof_values, std::vector< solution_symmetric_gradient_type< typename InputVector::value_type > > &symmetric_gradients) const
void get_function_gradients_from_local_dof_values(const InputVector &dof_values, std::vector< solution_gradient_type< typename InputVector::value_type > > &gradients) const
typename ProductType< Number, divergence_type >::type solution_divergence_type
void get_function_laplacians(const ReadVector< Number > &fe_function, std::vector< solution_laplacian_type< Number > > &laplacians) const
unsigned int first_vector_component
void get_function_gradients(const ReadVector< Number > &fe_function, std::vector< solution_gradient_type< Number > > &gradients) const
typename ProductType< Number, hessian_type >::type solution_hessian_type
void get_function_symmetric_gradients(const ReadVector< Number > &fe_function, std::vector< solution_symmetric_gradient_type< Number > > &symmetric_gradients) const
typename ProductType< Number, symmetric_gradient_type >::type solution_symmetric_gradient_type
void get_function_laplacians_from_local_dof_values(const InputVector &dof_values, std::vector< solution_laplacian_type< typename InputVector::value_type > > &laplacians) const
void get_function_hessians(const ReadVector< Number > &fe_function, std::vector< solution_hessian_type< Number > > &hessians) const
typename ProductType< Number, gradient_type >::type solution_gradient_type
typename ProductType< Number, value_type >::type solution_value_type
void get_function_third_derivatives(const ReadVector< Number > &fe_function, std::vector< solution_third_derivative_type< Number > > &third_derivatives) const
void get_function_values(const ReadVector< Number > &fe_function, std::vector< solution_value_type< Number > > &values) const
void get_function_curls_from_local_dof_values(const InputVector &dof_values, std::vector< solution_curl_type< typename InputVector::value_type > > &curls) const
typename ProductType< Number, curl_type >::type solution_curl_type
std::vector< ShapeFunctionData > shape_function_data
void get_function_curls(const ReadVector< Number > &fe_function, std::vector< solution_curl_type< Number > > &curls) const
void get_function_hessians_from_local_dof_values(const InputVector &dof_values, std::vector< solution_hessian_type< typename InputVector::value_type > > &hessians) const
typename ProductType< Number, value_type >::type solution_laplacian_type
void get_function_values_from_local_dof_values(const InputVector &dof_values, std::vector< solution_value_type< typename InputVector::value_type > > &values) const
void get_function_divergences_from_local_dof_values(const InputVector &dof_values, std::vector< solution_divergence_type< typename InputVector::value_type > > &divergences) const
void get_function_third_derivatives_from_local_dof_values(const InputVector &dof_values, std::vector< solution_third_derivative_type< typename InputVector::value_type > > &third_derivatives) const
void get_function_divergences(const ReadVector< Number > &fe_function, std::vector< solution_divergence_type< Number > > &divergences) const
unsigned int n_dofs_per_cell() const
unsigned int n_components() const
const ComponentMask & get_nonzero_components(const unsigned int i) const
bool is_primitive() const
std::pair< unsigned int, unsigned int > system_to_component_index(const unsigned int index) const
unsigned int n_nonzero_components(const unsigned int i) const
virtual size_type size() const =0
DEAL_II_HOST constexpr SymmetricTensor()=default
friend class Tensor
Definition tensor.h:865
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:518
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:519
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcMessage(std::string arg1)
@ update_hessians
Second derivatives of shape functions.
@ update_values
Shape function values.
@ update_3rd_derivatives
Third derivatives of shape functions.
@ update_gradients
Shape function gradients.
std::vector< unsigned int > make_shape_function_to_row_table(const FiniteElement< dim, spacedim > &fe)
Definition fe_values.cc:49
constexpr unsigned int invalid_unsigned_int
Definition types.h:232
Cache(const FEValuesBase< dim, spacedim > &fe_values)