deal.II version GIT relicensing-2238-gc05b561aad 2024-12-10 20:50: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\}}\)
Loading...
Searching...
No Matches
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
28
29
30namespace internal
31{
32 namespace
33 {
34 template <int dim, int spacedim>
35 inline std::vector<unsigned int>
37 {
38 std::vector<unsigned int> shape_function_to_row_table(
39 fe.n_dofs_per_cell() * fe.n_components(),
41 unsigned int row = 0;
42 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
43 {
44 // loop over all components that are nonzero for this particular
45 // shape function. if a component is zero then we leave the
46 // value in the table unchanged (at the invalid value)
47 // otherwise it is mapped to the next free entry
48 unsigned int nth_nonzero_component = 0;
49 for (unsigned int c = 0; c < fe.n_components(); ++c)
50 if (fe.get_nonzero_components(i)[c] == true)
51 {
52 shape_function_to_row_table[i * fe.n_components() + c] =
53 row + nth_nonzero_component;
54 ++nth_nonzero_component;
55 }
56 row += fe.n_nonzero_components(i);
57 }
58
59 return shape_function_to_row_table;
60 }
61 } // namespace
62} // namespace internal
63
64
65
66namespace FEValuesViews
67{
68 template <int dim, int spacedim>
70 const unsigned int component)
71 : fe_values(&fe_values)
72 , component(component)
73 , shape_function_data(this->fe_values->fe->n_dofs_per_cell())
74 {
75 const FiniteElement<dim, spacedim> &fe = *this->fe_values->fe;
77
78 // TODO: we'd like to use the fields with the same name as these
79 // variables from FEValuesBase, but they aren't initialized yet
80 // at the time we get here, so re-create it all
81 const std::vector<unsigned int> shape_function_to_row_table =
83
84 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
85 {
86 const bool is_primitive = fe.is_primitive() || fe.is_primitive(i);
87
88 if (is_primitive == true)
89 shape_function_data[i].is_nonzero_shape_function_component =
90 (component == fe.system_to_component_index(i).first);
91 else
92 shape_function_data[i].is_nonzero_shape_function_component =
93 (fe.get_nonzero_components(i)[component] == true);
94
95 if (shape_function_data[i].is_nonzero_shape_function_component == true)
96 shape_function_data[i].row_index =
97 shape_function_to_row_table[i * fe.n_components() + component];
98 else
100 }
101 }
102
103
104
105 template <int dim, int spacedim>
107 : fe_values(nullptr)
108 , component(numbers::invalid_unsigned_int)
109 {}
110
111
112
113 template <int dim, int spacedim>
115 const unsigned int first_vector_component)
116 : fe_values(&fe_values)
117 , first_vector_component(first_vector_component)
118 , shape_function_data(this->fe_values->fe->n_dofs_per_cell())
119 {
120 const FiniteElement<dim, spacedim> &fe = *this->fe_values->fe;
122
123 // TODO: we'd like to use the fields with the same name as these
124 // variables from FEValuesBase, but they aren't initialized yet
125 // at the time we get here, so re-create it all
126 const std::vector<unsigned int> shape_function_to_row_table =
128
129 for (unsigned int d = 0; d < spacedim; ++d)
130 {
131 const unsigned int component = first_vector_component + d;
132
133 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
134 {
135 const bool is_primitive = fe.is_primitive() || fe.is_primitive(i);
136
137 if (is_primitive == true)
138 shape_function_data[i].is_nonzero_shape_function_component[d] =
139 (component == fe.system_to_component_index(i).first);
140 else
141 shape_function_data[i].is_nonzero_shape_function_component[d] =
142 (fe.get_nonzero_components(i)[component] == true);
143
144 if (shape_function_data[i].is_nonzero_shape_function_component[d] ==
145 true)
146 shape_function_data[i].row_index[d] =
147 shape_function_to_row_table[i * fe.n_components() + component];
148 else
149 shape_function_data[i].row_index[d] =
151 }
152 }
153
154 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
155 {
156 unsigned int n_nonzero_components = 0;
157 for (unsigned int d = 0; d < spacedim; ++d)
158 if (shape_function_data[i].is_nonzero_shape_function_component[d] ==
159 true)
160 ++n_nonzero_components;
161
162 if (n_nonzero_components == 0)
163 shape_function_data[i].single_nonzero_component = -2;
164 else if (n_nonzero_components > 1)
165 shape_function_data[i].single_nonzero_component = -1;
166 else
167 {
168 for (unsigned int d = 0; d < spacedim; ++d)
170 .is_nonzero_shape_function_component[d] == true)
171 {
172 shape_function_data[i].single_nonzero_component =
173 shape_function_data[i].row_index[d];
174 shape_function_data[i].single_nonzero_component_index = d;
175 break;
176 }
177 }
178 }
179 }
180
181
182
183 template <int dim, int spacedim>
185 : fe_values(nullptr)
186 , first_vector_component(numbers::invalid_unsigned_int)
187 {}
188
189
190
191 template <int dim, int spacedim>
193 const FEValuesBase<dim, spacedim> &fe_values,
194 const unsigned int first_tensor_component)
195 : fe_values(&fe_values)
196 , first_tensor_component(first_tensor_component)
197 , shape_function_data(this->fe_values->fe->n_dofs_per_cell())
198 {
199 const FiniteElement<dim, spacedim> &fe = *this->fe_values->fe;
200 Assert(first_tensor_component + (dim * dim + dim) / 2 - 1 <
201 fe.n_components(),
203 first_tensor_component +
205 0,
206 fe.n_components()));
207 // TODO: we'd like to use the fields with the same name as these
208 // variables from FEValuesBase, but they aren't initialized yet
209 // at the time we get here, so re-create it all
210 const std::vector<unsigned int> shape_function_to_row_table =
212
213 for (unsigned int d = 0;
214 d < ::SymmetricTensor<2, dim>::n_independent_components;
215 ++d)
216 {
217 const unsigned int component = first_tensor_component + d;
218
219 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
220 {
221 const bool is_primitive = fe.is_primitive() || fe.is_primitive(i);
222
223 if (is_primitive == true)
224 shape_function_data[i].is_nonzero_shape_function_component[d] =
225 (component == fe.system_to_component_index(i).first);
226 else
227 shape_function_data[i].is_nonzero_shape_function_component[d] =
228 (fe.get_nonzero_components(i)[component] == true);
229
230 if (shape_function_data[i].is_nonzero_shape_function_component[d] ==
231 true)
232 shape_function_data[i].row_index[d] =
233 shape_function_to_row_table[i * fe.n_components() + component];
234 else
235 shape_function_data[i].row_index[d] =
237 }
238 }
239
240 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
241 {
242 unsigned int n_nonzero_components = 0;
243 for (unsigned int d = 0;
244 d < ::SymmetricTensor<2, dim>::n_independent_components;
245 ++d)
246 if (shape_function_data[i].is_nonzero_shape_function_component[d] ==
247 true)
248 ++n_nonzero_components;
249
250 if (n_nonzero_components == 0)
251 shape_function_data[i].single_nonzero_component = -2;
252 else if (n_nonzero_components > 1)
253 shape_function_data[i].single_nonzero_component = -1;
254 else
255 {
256 for (unsigned int d = 0;
257 d < ::SymmetricTensor<2, dim>::n_independent_components;
258 ++d)
259 if (shape_function_data[i]
260 .is_nonzero_shape_function_component[d] == true)
261 {
262 shape_function_data[i].single_nonzero_component =
263 shape_function_data[i].row_index[d];
264 shape_function_data[i].single_nonzero_component_index = d;
265 break;
266 }
267 }
268 }
269 }
270
271
272
273 template <int dim, int spacedim>
275 : fe_values(nullptr)
276 , first_tensor_component(numbers::invalid_unsigned_int)
277 {}
278
279
280
281 template <int dim, int spacedim>
283 const unsigned int first_tensor_component)
284 : fe_values(&fe_values)
285 , first_tensor_component(first_tensor_component)
286 , shape_function_data(this->fe_values->fe->n_dofs_per_cell())
287 {
288 const FiniteElement<dim, spacedim> &fe = *this->fe_values->fe;
289 AssertIndexRange(first_tensor_component + dim * dim - 1, fe.n_components());
290 // TODO: we'd like to use the fields with the same name as these
291 // variables from FEValuesBase, but they aren't initialized yet
292 // at the time we get here, so re-create it all
293 const std::vector<unsigned int> shape_function_to_row_table =
295
296 for (unsigned int d = 0; d < dim * dim; ++d)
297 {
298 const unsigned int component = first_tensor_component + d;
299
300 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
301 {
302 const bool is_primitive = fe.is_primitive() || fe.is_primitive(i);
303
304 if (is_primitive == true)
305 shape_function_data[i].is_nonzero_shape_function_component[d] =
306 (component == fe.system_to_component_index(i).first);
307 else
308 shape_function_data[i].is_nonzero_shape_function_component[d] =
309 (fe.get_nonzero_components(i)[component] == true);
310
311 if (shape_function_data[i].is_nonzero_shape_function_component[d] ==
312 true)
313 shape_function_data[i].row_index[d] =
314 shape_function_to_row_table[i * fe.n_components() + component];
315 else
316 shape_function_data[i].row_index[d] =
318 }
319 }
320
321 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
322 {
323 unsigned int n_nonzero_components = 0;
324 for (unsigned int d = 0; d < dim * dim; ++d)
325 if (shape_function_data[i].is_nonzero_shape_function_component[d] ==
326 true)
327 ++n_nonzero_components;
328
329 if (n_nonzero_components == 0)
330 shape_function_data[i].single_nonzero_component = -2;
331 else if (n_nonzero_components > 1)
332 shape_function_data[i].single_nonzero_component = -1;
333 else
334 {
335 for (unsigned int d = 0; d < dim * dim; ++d)
336 if (shape_function_data[i]
337 .is_nonzero_shape_function_component[d] == true)
338 {
339 shape_function_data[i].single_nonzero_component =
340 shape_function_data[i].row_index[d];
341 shape_function_data[i].single_nonzero_component_index = d;
342 break;
343 }
344 }
345 }
346 }
347
348
349
350 template <int dim, int spacedim>
352 : fe_values(nullptr)
353 , first_tensor_component(numbers::invalid_unsigned_int)
354 {}
355
356
357
358 template <int dim, int spacedim>
359 template <typename Number>
360 void
362 const ReadVector<Number> &fe_function,
363 std::vector<solution_value_type<Number>> &values) const
364 {
365 Assert(fe_values->update_flags & update_values,
367 "update_values")));
368 Assert(fe_values->present_cell.is_initialized(),
370 AssertDimension(values.size(), fe_values->n_quadrature_points);
371
372 // get function values of dofs on this cell and call internal worker
373 // function
374 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
375 fe_values->present_cell.get_interpolated_dof_values(fe_function,
376 dof_values);
377 internal::do_function_values<dim, spacedim>(
378 make_const_array_view(dof_values),
379 fe_values->finite_element_output.shape_values,
380 shape_function_data,
381 values);
382 }
383
384
385
386 template <int dim, int spacedim>
387 template <class InputVector>
388 void
390 const InputVector &dof_values,
392 const
393 {
394 Assert(fe_values->update_flags & update_values,
396 "update_values")));
397 Assert(fe_values->present_cell.is_initialized(),
399 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
400 AssertDimension(values.size(), fe_values->n_quadrature_points);
401
402 internal::do_function_values<dim, spacedim>(
403 make_const_array_view(dof_values),
404 fe_values->finite_element_output.shape_values,
405 shape_function_data,
406 values);
407 }
408
409
410
411 template <int dim, int spacedim>
412 template <typename Number>
413 void
415 const ReadVector<Number> &fe_function,
416 std::vector<solution_gradient_type<Number>> &gradients) const
417 {
418 Assert(fe_values->update_flags & update_gradients,
420 "update_gradients")));
421 Assert(fe_values->present_cell.is_initialized(),
423 AssertDimension(fe_function.size(),
424 fe_values->present_cell.n_dofs_for_dof_handler());
425 AssertDimension(gradients.size(), fe_values->n_quadrature_points);
426
427 // get function values of dofs on this cell
428 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
429 fe_values->present_cell.get_interpolated_dof_values(fe_function,
430 dof_values);
431 internal::do_function_derivatives<1, dim, spacedim>(
432 make_const_array_view(dof_values),
433 fe_values->finite_element_output.shape_gradients,
434 shape_function_data,
435 gradients);
436 }
437
438
439
440 template <int dim, int spacedim>
441 template <typename InputVector>
442 void
444 const InputVector &dof_values,
446 &gradients) const
447 {
448 Assert(fe_values->update_flags & update_gradients,
450 "update_gradients")));
451 Assert(fe_values->present_cell.is_initialized(),
453 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
454 AssertDimension(gradients.size(), fe_values->n_quadrature_points);
455
456 internal::do_function_derivatives<1, dim, spacedim>(
457 make_const_array_view(dof_values),
458 fe_values->finite_element_output.shape_gradients,
459 shape_function_data,
460 gradients);
461 }
462
463
464
465 template <int dim, int spacedim>
466 template <typename Number>
467 void
469 const ReadVector<Number> &fe_function,
470 std::vector<solution_hessian_type<Number>> &hessians) const
471 {
472 Assert(fe_values->update_flags & update_hessians,
474 "update_hessians")));
475 Assert(fe_values->present_cell.is_initialized(),
477 AssertDimension(fe_function.size(),
478 fe_values->present_cell.n_dofs_for_dof_handler());
479 AssertDimension(hessians.size(), fe_values->n_quadrature_points);
480
481 // get function values of dofs on this cell
482 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
483 fe_values->present_cell.get_interpolated_dof_values(fe_function,
484 dof_values);
485 internal::do_function_derivatives<2, dim, spacedim>(
486 make_const_array_view(dof_values),
487 fe_values->finite_element_output.shape_hessians,
488 shape_function_data,
489 hessians);
490 }
491
492
493
494 template <int dim, int spacedim>
495 template <class InputVector>
496 void
498 const InputVector &dof_values,
500 &hessians) const
501 {
502 Assert(fe_values->update_flags & update_hessians,
504 "update_hessians")));
505 Assert(fe_values->present_cell.is_initialized(),
507 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
508 AssertDimension(hessians.size(), fe_values->n_quadrature_points);
509
510 internal::do_function_derivatives<2, dim, spacedim>(
511 make_const_array_view(dof_values),
512 fe_values->finite_element_output.shape_hessians,
513 shape_function_data,
514 hessians);
515 }
516
517
518
519 template <int dim, int spacedim>
520 template <typename Number>
521 void
523 const ReadVector<Number> &fe_function,
524 std::vector<solution_laplacian_type<Number>> &laplacians) const
525 {
526 Assert(fe_values->update_flags & update_hessians,
528 "update_hessians")));
529 Assert(fe_values->present_cell.is_initialized(),
531 AssertDimension(fe_function.size(),
532 fe_values->present_cell.n_dofs_for_dof_handler());
533 AssertDimension(laplacians.size(), fe_values->n_quadrature_points);
534
535 // get function values of dofs on this cell
536 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
537 fe_values->present_cell.get_interpolated_dof_values(fe_function,
538 dof_values);
539 internal::do_function_laplacians<dim, spacedim>(
540 make_const_array_view(dof_values),
541 fe_values->finite_element_output.shape_hessians,
542 shape_function_data,
543 laplacians);
544 }
545
546
547
548 template <int dim, int spacedim>
549 template <class InputVector>
550 void
552 const InputVector &dof_values,
554 &laplacians) const
555 {
556 Assert(fe_values->update_flags & update_hessians,
558 "update_hessians")));
559 Assert(fe_values->present_cell.is_initialized(),
561 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
562 AssertDimension(laplacians.size(), fe_values->n_quadrature_points);
563
564 internal::do_function_laplacians<dim, spacedim>(
565 make_const_array_view(dof_values),
566 fe_values->finite_element_output.shape_hessians,
567 shape_function_data,
568 laplacians);
569 }
570
571
572
573 template <int dim, int spacedim>
574 template <typename Number>
575 void
577 const ReadVector<Number> &fe_function,
578 std::vector<solution_third_derivative_type<Number>> &third_derivatives)
579 const
580 {
581 Assert(fe_values->update_flags & update_3rd_derivatives,
583 "update_3rd_derivatives")));
584 Assert(fe_values->present_cell.is_initialized(),
586 AssertDimension(fe_function.size(),
587 fe_values->present_cell.n_dofs_for_dof_handler());
588 AssertDimension(third_derivatives.size(), fe_values->n_quadrature_points);
589
590 // get function values of dofs on this cell
591 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
592 fe_values->present_cell.get_interpolated_dof_values(fe_function,
593 dof_values);
594 internal::do_function_derivatives<3, dim, spacedim>(
595 make_const_array_view(dof_values),
596 fe_values->finite_element_output.shape_3rd_derivatives,
597 shape_function_data,
598 third_derivatives);
599 }
600
601
602
603 template <int dim, int spacedim>
604 template <class InputVector>
605 void
607 const InputVector &dof_values,
608 std::vector<
610 &third_derivatives) const
611 {
612 Assert(fe_values->update_flags & update_3rd_derivatives,
614 "update_3rd_derivatives")));
615 Assert(fe_values->present_cell.is_initialized(),
617 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
618 AssertDimension(third_derivatives.size(), fe_values->n_quadrature_points);
619
620 internal::do_function_derivatives<3, dim, spacedim>(
621 make_const_array_view(dof_values),
622 fe_values->finite_element_output.shape_3rd_derivatives,
623 shape_function_data,
624 third_derivatives);
625 }
626
627
628
629 template <int dim, int spacedim>
630 template <typename Number>
631 void
633 const ReadVector<Number> &fe_function,
634 std::vector<solution_value_type<Number>> &values) const
635 {
636 Assert(fe_values->update_flags & update_values,
638 "update_values")));
639 Assert(fe_values->present_cell.is_initialized(),
641 AssertDimension(values.size(), fe_values->n_quadrature_points);
642
643 // get function values of dofs on this cell
644 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
645 fe_values->present_cell.get_interpolated_dof_values(fe_function,
646 dof_values);
647 internal::do_function_values<dim, spacedim>(
648 make_const_array_view(dof_values),
649 fe_values->finite_element_output.shape_values,
650 shape_function_data,
651 values);
652 }
653
654
655
656 template <int dim, int spacedim>
657 template <class InputVector>
658 void
660 const InputVector &dof_values,
662 const
663 {
664 Assert(fe_values->update_flags & update_values,
666 "update_values")));
667 Assert(fe_values->present_cell.is_initialized(),
669 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
670 AssertDimension(values.size(), fe_values->n_quadrature_points);
671
672 internal::do_function_values<dim, spacedim>(
673 make_const_array_view(dof_values),
674 fe_values->finite_element_output.shape_values,
675 shape_function_data,
676 values);
677 }
678
679
680
681 template <int dim, int spacedim>
682 template <typename Number>
683 void
685 const ReadVector<Number> &fe_function,
686 std::vector<solution_gradient_type<Number>> &gradients) const
687 {
688 Assert(fe_values->update_flags & update_gradients,
690 "update_gradients")));
691 Assert(fe_values->present_cell.is_initialized(),
693 AssertDimension(fe_function.size(),
694 fe_values->present_cell.n_dofs_for_dof_handler());
695 AssertDimension(gradients.size(), fe_values->n_quadrature_points);
696
697 // get function values of dofs on this cell
698 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
699 fe_values->present_cell.get_interpolated_dof_values(fe_function,
700 dof_values);
701 internal::do_function_derivatives<1, dim, spacedim>(
702 make_const_array_view(dof_values),
703 fe_values->finite_element_output.shape_gradients,
704 shape_function_data,
705 gradients);
706 }
707
708
709
710 template <int dim, int spacedim>
711 template <typename InputVector>
712 void
714 const InputVector &dof_values,
716 &gradients) const
717 {
718 Assert(fe_values->update_flags & update_gradients,
720 "update_gradients")));
721 Assert(fe_values->present_cell.is_initialized(),
723 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
724 AssertDimension(gradients.size(), fe_values->n_quadrature_points);
725
726 internal::do_function_derivatives<1, dim, spacedim>(
727 make_const_array_view(dof_values),
728 fe_values->finite_element_output.shape_gradients,
729 shape_function_data,
730 gradients);
731 }
732
733
734
735 template <int dim, int spacedim>
736 template <typename Number>
737 void
739 const ReadVector<Number> &fe_function,
740 std::vector<solution_symmetric_gradient_type<Number>> &symmetric_gradients)
741 const
742 {
743 Assert(fe_values->update_flags & update_gradients,
745 "update_gradients")));
746 Assert(fe_values->present_cell.is_initialized(),
748 AssertDimension(fe_function.size(),
749 fe_values->present_cell.n_dofs_for_dof_handler());
750 AssertDimension(symmetric_gradients.size(), fe_values->n_quadrature_points);
751
752 // get function values of dofs on this cell
753 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
754 fe_values->present_cell.get_interpolated_dof_values(fe_function,
755 dof_values);
756 internal::do_function_symmetric_gradients<dim, spacedim>(
757 make_const_array_view(dof_values),
758 fe_values->finite_element_output.shape_gradients,
759 shape_function_data,
760 symmetric_gradients);
761 }
762
763
764
765 template <int dim, int spacedim>
766 template <class InputVector>
767 void
769 const InputVector &dof_values,
770 std::vector<
772 &symmetric_gradients) const
773 {
774 Assert(fe_values->update_flags & update_gradients,
776 "update_gradients")));
777 Assert(fe_values->present_cell.is_initialized(),
779 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
780 AssertDimension(symmetric_gradients.size(), fe_values->n_quadrature_points);
781
782 internal::do_function_symmetric_gradients<dim, spacedim>(
783 make_const_array_view(dof_values),
784 fe_values->finite_element_output.shape_gradients,
785 shape_function_data,
786 symmetric_gradients);
787 }
788
789
790
791 template <int dim, int spacedim>
792 template <typename Number>
793 void
795 const ReadVector<Number> &fe_function,
796 std::vector<solution_divergence_type<Number>> &divergences) const
797 {
798 Assert(fe_values->update_flags & update_gradients,
800 "update_gradients")));
801 Assert(fe_values->present_cell.is_initialized(),
803 AssertDimension(fe_function.size(),
804 fe_values->present_cell.n_dofs_for_dof_handler());
805 AssertDimension(divergences.size(), fe_values->n_quadrature_points);
806
807 // get function values of dofs
808 // on this cell
809 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
810 fe_values->present_cell.get_interpolated_dof_values(fe_function,
811 dof_values);
812 internal::do_function_divergences<dim, spacedim>(
813 make_const_array_view(dof_values),
814 fe_values->finite_element_output.shape_gradients,
815 shape_function_data,
816 divergences);
817 }
818
819
820
821 template <int dim, int spacedim>
822 template <class InputVector>
823 void
825 const InputVector &dof_values,
827 &divergences) const
828 {
829 Assert(fe_values->update_flags & update_gradients,
831 "update_gradients")));
832 Assert(fe_values->present_cell.is_initialized(),
834 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
835 AssertDimension(divergences.size(), fe_values->n_quadrature_points);
836
837 internal::do_function_divergences<dim, spacedim>(
838 make_const_array_view(dof_values),
839 fe_values->finite_element_output.shape_gradients,
840 shape_function_data,
841 divergences);
842 }
843
844
845
846 template <int dim, int spacedim>
847 template <typename Number>
848 void
850 const ReadVector<Number> &fe_function,
851 std::vector<solution_curl_type<Number>> &curls) const
852 {
853 Assert(fe_values->update_flags & update_gradients,
855 "update_gradients")));
856 Assert(fe_values->present_cell.is_initialized(),
857 ExcMessage("FEValues object is not reinited to any cell"));
858 AssertDimension(fe_function.size(),
859 fe_values->present_cell.n_dofs_for_dof_handler());
860 AssertDimension(curls.size(), fe_values->n_quadrature_points);
861
862 // get function values of dofs on this cell
863 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
864 fe_values->present_cell.get_interpolated_dof_values(fe_function,
865 dof_values);
866 internal::do_function_curls<dim, spacedim>(
867 make_const_array_view(dof_values),
868 fe_values->finite_element_output.shape_gradients,
869 shape_function_data,
870 curls);
871 }
872
873
874
875 template <int dim, int spacedim>
876 template <class InputVector>
877 void
879 const InputVector &dof_values,
881 const
882 {
883 Assert(fe_values->update_flags & update_gradients,
885 "update_gradients")));
886 Assert(fe_values->present_cell.is_initialized(),
887 ExcMessage("FEValues object is not reinited to any cell"));
888 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
889 AssertDimension(curls.size(), fe_values->n_quadrature_points);
890
891 internal::do_function_curls<dim, spacedim>(
892 make_const_array_view(dof_values),
893 fe_values->finite_element_output.shape_gradients,
894 shape_function_data,
895 curls);
896 }
897
898
899
900 template <int dim, int spacedim>
901 template <typename Number>
902 void
904 const ReadVector<Number> &fe_function,
905 std::vector<solution_hessian_type<Number>> &hessians) const
906 {
907 Assert(fe_values->update_flags & update_hessians,
909 "update_hessians")));
910 Assert(fe_values->present_cell.is_initialized(),
912 AssertDimension(fe_function.size(),
913 fe_values->present_cell.n_dofs_for_dof_handler());
914 AssertDimension(hessians.size(), fe_values->n_quadrature_points);
915
916 // get function values of dofs on this cell
917 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
918 fe_values->present_cell.get_interpolated_dof_values(fe_function,
919 dof_values);
920 internal::do_function_derivatives<2, dim, spacedim>(
921 make_const_array_view(dof_values),
922 fe_values->finite_element_output.shape_hessians,
923 shape_function_data,
924 hessians);
925 }
926
927
928
929 template <int dim, int spacedim>
930 template <class InputVector>
931 void
933 const InputVector &dof_values,
935 &hessians) const
936 {
937 Assert(fe_values->update_flags & update_hessians,
939 "update_hessians")));
940 Assert(fe_values->present_cell.is_initialized(),
942 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
943 AssertDimension(hessians.size(), fe_values->n_quadrature_points);
944
945 internal::do_function_derivatives<2, dim, spacedim>(
946 make_const_array_view(dof_values),
947 fe_values->finite_element_output.shape_hessians,
948 shape_function_data,
949 hessians);
950 }
951
952
953
954 template <int dim, int spacedim>
955 template <typename Number>
956 void
958 const ReadVector<Number> &fe_function,
959 std::vector<solution_value_type<Number>> &laplacians) const
960 {
961 Assert(fe_values->update_flags & update_hessians,
963 "update_hessians")));
964 Assert(laplacians.size() == fe_values->n_quadrature_points,
965 ExcDimensionMismatch(laplacians.size(),
966 fe_values->n_quadrature_points));
967 Assert(fe_values->present_cell.is_initialized(),
969 Assert(
970 fe_function.size() == fe_values->present_cell.n_dofs_for_dof_handler(),
971 ExcDimensionMismatch(fe_function.size(),
972 fe_values->present_cell.n_dofs_for_dof_handler()));
973 AssertDimension(laplacians.size(), fe_values->n_quadrature_points);
974
975 // get function values of dofs on this cell
976 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
977 fe_values->present_cell.get_interpolated_dof_values(fe_function,
978 dof_values);
979 internal::do_function_laplacians<dim, spacedim>(
980 make_const_array_view(dof_values),
981 fe_values->finite_element_output.shape_hessians,
982 shape_function_data,
983 laplacians);
984 }
985
986
987
988 template <int dim, int spacedim>
989 template <class InputVector>
990 void
992 const InputVector &dof_values,
994 &laplacians) const
995 {
996 Assert(fe_values->update_flags & update_hessians,
998 "update_hessians")));
999 Assert(laplacians.size() == fe_values->n_quadrature_points,
1000 ExcDimensionMismatch(laplacians.size(),
1001 fe_values->n_quadrature_points));
1002 Assert(fe_values->present_cell.is_initialized(),
1004 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
1005 AssertDimension(laplacians.size(), fe_values->n_quadrature_points);
1006
1007 internal::do_function_laplacians<dim, spacedim>(
1008 make_const_array_view(dof_values),
1009 fe_values->finite_element_output.shape_hessians,
1010 shape_function_data,
1011 laplacians);
1012 }
1013
1014
1015
1016 template <int dim, int spacedim>
1017 template <typename Number>
1018 void
1020 const ReadVector<Number> &fe_function,
1021 std::vector<solution_third_derivative_type<Number>> &third_derivatives)
1022 const
1023 {
1024 Assert(fe_values->update_flags & update_3rd_derivatives,
1026 "update_3rd_derivatives")));
1027 Assert(fe_values->present_cell.is_initialized(),
1029 AssertDimension(fe_function.size(),
1030 fe_values->present_cell.n_dofs_for_dof_handler());
1031 AssertDimension(third_derivatives.size(), fe_values->n_quadrature_points);
1032
1033 // get function values of dofs on this cell
1034 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
1035 fe_values->present_cell.get_interpolated_dof_values(fe_function,
1036 dof_values);
1037 internal::do_function_derivatives<3, dim, spacedim>(
1038 make_const_array_view(dof_values),
1039 fe_values->finite_element_output.shape_3rd_derivatives,
1040 shape_function_data,
1041 third_derivatives);
1042 }
1043
1044
1045
1046 template <int dim, int spacedim>
1047 template <class InputVector>
1048 void
1050 const InputVector &dof_values,
1051 std::vector<
1053 &third_derivatives) const
1054 {
1055 Assert(fe_values->update_flags & update_3rd_derivatives,
1057 "update_3rd_derivatives")));
1058 Assert(fe_values->present_cell.is_initialized(),
1060 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
1061 AssertDimension(third_derivatives.size(), fe_values->n_quadrature_points);
1062
1063 internal::do_function_derivatives<3, dim, spacedim>(
1064 make_const_array_view(dof_values),
1065 fe_values->finite_element_output.shape_3rd_derivatives,
1066 shape_function_data,
1067 third_derivatives);
1068 }
1069
1070
1071
1072 template <int dim, int spacedim>
1073 template <typename Number>
1074 void
1076 const ReadVector<Number> &fe_function,
1077 std::vector<solution_value_type<Number>> &values) const
1078 {
1079 Assert(fe_values->update_flags & update_values,
1081 "update_values")));
1082 Assert(fe_values->present_cell.is_initialized(),
1084 AssertDimension(values.size(), fe_values->n_quadrature_points);
1085
1086 // get function values of dofs on this cell
1087 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
1088 fe_values->present_cell.get_interpolated_dof_values(fe_function,
1089 dof_values);
1090 internal::do_function_values<dim, spacedim>(
1091 make_const_array_view(dof_values),
1092 fe_values->finite_element_output.shape_values,
1093 shape_function_data,
1094 values);
1095 }
1096
1097
1098
1099 template <int dim, int spacedim>
1100 template <class InputVector>
1101 void
1103 const InputVector &dof_values,
1105 const
1106 {
1107 Assert(fe_values->update_flags & update_values,
1109 "update_values")));
1110 Assert(fe_values->present_cell.is_initialized(),
1112 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
1113 AssertDimension(values.size(), fe_values->n_quadrature_points);
1114
1115 internal::do_function_values<dim, spacedim>(
1116 make_const_array_view(dof_values),
1117 fe_values->finite_element_output.shape_values,
1118 shape_function_data,
1119 values);
1120 }
1121
1122
1123
1124 template <int dim, int spacedim>
1125 template <typename Number>
1126 void
1128 const ReadVector<Number> &fe_function,
1129 std::vector<solution_divergence_type<Number>> &divergences) const
1130 {
1131 Assert(fe_values->update_flags & update_gradients,
1133 "update_gradients")));
1134 Assert(fe_values->present_cell.is_initialized(),
1136 AssertDimension(fe_function.size(),
1137 fe_values->present_cell.n_dofs_for_dof_handler());
1138 AssertDimension(divergences.size(), fe_values->n_quadrature_points);
1139
1140 // get function values of dofs
1141 // on this cell
1142 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
1143 fe_values->present_cell.get_interpolated_dof_values(fe_function,
1144 dof_values);
1145 internal::do_function_divergences<dim, spacedim>(
1146 make_const_array_view(dof_values),
1147 fe_values->finite_element_output.shape_gradients,
1148 shape_function_data,
1149 divergences);
1150 }
1151
1152
1153
1154 template <int dim, int spacedim>
1155 template <class InputVector>
1156 void
1159 const InputVector &dof_values,
1161 &divergences) const
1162 {
1163 Assert(fe_values->update_flags & update_gradients,
1165 "update_gradients")));
1166 Assert(fe_values->present_cell.is_initialized(),
1168 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
1169 AssertDimension(divergences.size(), fe_values->n_quadrature_points);
1170
1171 internal::do_function_divergences<dim, spacedim>(
1172 make_const_array_view(dof_values),
1173 fe_values->finite_element_output.shape_gradients,
1174 shape_function_data,
1175 divergences);
1176 }
1177
1178
1179
1180 template <int dim, int spacedim>
1181 template <typename Number>
1182 void
1184 const ReadVector<Number> &fe_function,
1185 std::vector<solution_value_type<Number>> &values) const
1186 {
1187 Assert(fe_values->update_flags & update_values,
1189 "update_values")));
1190 Assert(fe_values->present_cell.is_initialized(),
1192 AssertDimension(values.size(), fe_values->n_quadrature_points);
1193
1194 // get function values of dofs on this cell
1195 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
1196 fe_values->present_cell.get_interpolated_dof_values(fe_function,
1197 dof_values);
1198 internal::do_function_values<dim, spacedim>(
1199 make_const_array_view(dof_values),
1200 fe_values->finite_element_output.shape_values,
1201 shape_function_data,
1202 values);
1203 }
1204
1205
1206
1207 template <int dim, int spacedim>
1208 template <class InputVector>
1209 void
1211 const InputVector &dof_values,
1213 const
1214 {
1215 Assert(fe_values->update_flags & update_values,
1217 "update_values")));
1218 Assert(fe_values->present_cell.is_initialized(),
1220 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
1221 AssertDimension(values.size(), fe_values->n_quadrature_points);
1222
1223 internal::do_function_values<dim, spacedim>(
1224 make_const_array_view(dof_values),
1225 fe_values->finite_element_output.shape_values,
1226 shape_function_data,
1227 values);
1228 }
1229
1230
1231
1232 template <int dim, int spacedim>
1233 template <typename Number>
1234 void
1236 const ReadVector<Number> &fe_function,
1237 std::vector<solution_divergence_type<Number>> &divergences) const
1238 {
1239 Assert(fe_values->update_flags & update_gradients,
1241 "update_gradients")));
1242 Assert(fe_values->present_cell.is_initialized(),
1244 AssertDimension(fe_function.size(),
1245 fe_values->present_cell.n_dofs_for_dof_handler());
1246 AssertDimension(divergences.size(), fe_values->n_quadrature_points);
1247
1248 // get function values of dofs
1249 // on this cell
1250 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
1251 fe_values->present_cell.get_interpolated_dof_values(fe_function,
1252 dof_values);
1253 internal::do_function_divergences<dim, spacedim>(
1254 make_const_array_view(dof_values),
1255 fe_values->finite_element_output.shape_gradients,
1256 shape_function_data,
1257 divergences);
1258 }
1259
1260
1261
1262 template <int dim, int spacedim>
1263 template <class InputVector>
1264 void
1266 const InputVector &dof_values,
1268 &divergences) const
1269 {
1270 Assert(fe_values->update_flags & update_gradients,
1272 "update_gradients")));
1273 Assert(fe_values->present_cell.is_initialized(),
1275 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
1276 AssertDimension(divergences.size(), fe_values->n_quadrature_points);
1277
1278 internal::do_function_divergences<dim, spacedim>(
1279 make_const_array_view(dof_values),
1280 fe_values->finite_element_output.shape_gradients,
1281 shape_function_data,
1282 divergences);
1283 }
1284
1285
1286
1287 template <int dim, int spacedim>
1288 template <typename Number>
1289 void
1291 const ReadVector<Number> &fe_function,
1292 std::vector<solution_gradient_type<Number>> &gradients) const
1293 {
1294 Assert(fe_values->update_flags & update_gradients,
1296 "update_gradients")));
1297 Assert(fe_values->present_cell.is_initialized(),
1299 AssertDimension(fe_function.size(),
1300 fe_values->present_cell.n_dofs_for_dof_handler());
1301 AssertDimension(gradients.size(), fe_values->n_quadrature_points);
1302
1303 // get function values of dofs
1304 // on this cell
1305 ::Vector<Number> dof_values(fe_values->dofs_per_cell);
1306 fe_values->present_cell.get_interpolated_dof_values(fe_function,
1307 dof_values);
1308 internal::do_function_gradients<dim, spacedim>(
1309 make_const_array_view(dof_values),
1310 fe_values->finite_element_output.shape_gradients,
1311 shape_function_data,
1312 gradients);
1313 }
1314
1315
1316
1317 template <int dim, int spacedim>
1318 template <class InputVector>
1319 void
1321 const InputVector &dof_values,
1323 &gradients) const
1324 {
1325 Assert(fe_values->update_flags & update_gradients,
1327 "update_gradients")));
1328 Assert(fe_values->present_cell.is_initialized(),
1330 AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
1331 AssertDimension(gradients.size(), fe_values->n_quadrature_points);
1332
1333 internal::do_function_gradients<dim, spacedim>(
1334 make_const_array_view(dof_values),
1335 fe_values->finite_element_output.shape_gradients,
1336 shape_function_data,
1337 gradients);
1338 }
1339} // namespace FEValuesViews
1340
1341
1342namespace internal
1343{
1344 namespace FEValuesViews
1345 {
1346 template <int dim, int spacedim>
1348 {
1349 const FiniteElement<dim, spacedim> &fe = fe_values.get_fe();
1350
1351 const unsigned int n_scalars = fe.n_components();
1352 scalars.resize(n_scalars);
1353
1354 // compute number of vectors that we can fit into this finite element.
1355 // note that this is based on the dimensionality 'dim' of the manifold,
1356 // not 'spacedim' of the output vector
1357 const unsigned int n_vectors =
1360 1 :
1361 0);
1362 vectors.resize(n_vectors);
1363
1364 // compute number of symmetric tensors in the same way as above
1365 const unsigned int n_symmetric_second_order_tensors =
1366 (fe.n_components() >=
1368 fe.n_components() -
1370 0);
1371 symmetric_second_order_tensors.resize(n_symmetric_second_order_tensors);
1372
1373 // compute number of symmetric tensors in the same way as above
1374 const unsigned int n_second_order_tensors =
1377 1 :
1378 0);
1379 second_order_tensors.resize(n_second_order_tensors);
1380 }
1381 } // namespace FEValuesViews
1382} // namespace internal
1383
1384/*------------------------------- Explicit Instantiations -------------*/
1385
1386#include "fe_values_views.inst"
1387
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:498
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:499
#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
static const unsigned int invalid_unsigned_int
Definition types.h:220
Cache(const FEValuesBase< dim, spacedim > &fe_values)