Loading [MathJax]/extensions/TeX/newcommand.js
 Reference documentation for deal.II version 9.3.3
\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
mapping_cartesian.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2001 - 2021 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
21#include <deal.II/base/tensor.h>
22
24
27
28#include <deal.II/grid/tria.h>
30
32
33#include <algorithm>
34#include <cmath>
35#include <memory>
36
37
39
40
41
42template <int dim, int spacedim>
44 const Quadrature<dim> &q)
45 : cell_extents(numbers::signaling_nan<Tensor<1, dim>>())
46 , volume_element(numbers::signaling_nan<double>())
47 , quadrature_points(q.get_points())
48{}
49
50
51
52template <int dim, int spacedim>
53std::size_t
55{
60}
61
62
63
64template <int dim, int spacedim>
65bool
67{
68 return true;
69}
70
71
72
73template <int dim, int spacedim>
74bool
76 const ReferenceCell &reference_cell) const
77{
78 Assert(dim == reference_cell.get_dimension(),
79 ExcMessage("The dimension of your mapping (" +
81 ") and the reference cell cell_type (" +
82 Utilities::to_string(reference_cell.get_dimension()) +
83 " ) do not agree."));
84
85 return reference_cell.is_hyper_cube();
86}
87
88
89
90template <int dim, int spacedim>
93 const UpdateFlags in) const
94{
95 // this mapping is pretty simple in that it can basically compute
96 // every piece of information wanted by FEValues without requiring
97 // computing any other quantities. boundary forms are one exception
98 // since they can be computed from the normal vectors without much
99 // further ado
100 UpdateFlags out = in;
101 if (out & update_boundary_forms)
103
104 return out;
105}
106
107
108
109template <int dim, int spacedim>
110std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
112 const Quadrature<dim> &q) const
113{
114 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
115 std::make_unique<InternalData>(q);
116 auto &data = dynamic_cast<InternalData &>(*data_ptr);
117
118 // store the flags in the internal data object so we can access them
119 // in fill_fe_*_values(). use the transitive hull of the required
120 // flags
121 data.update_each = requires_update_flags(update_flags);
122
123 return data_ptr;
124}
125
126
127
128template <int dim, int spacedim>
129std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
131 const UpdateFlags update_flags,
132 const hp::QCollection<dim - 1> &quadrature) const
133{
134 AssertDimension(quadrature.size(), 1);
135
136 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
137 std::make_unique<InternalData>(QProjector<dim>::project_to_all_faces(
138 ReferenceCells::get_hypercube<dim>(), quadrature[0]));
139 auto &data = dynamic_cast<InternalData &>(*data_ptr);
140
141 // verify that we have computed the transitive hull of the required
142 // flags and that FEValues has faithfully passed them on to us
143 Assert(update_flags == requires_update_flags(update_flags),
145
146 // store the flags in the internal data object so we can access them
147 // in fill_fe_*_values()
148 data.update_each = update_flags;
149
150 return data_ptr;
151}
152
153
154
155template <int dim, int spacedim>
156std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
158 const UpdateFlags update_flags,
159 const Quadrature<dim - 1> &quadrature) const
160{
161 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
162 std::make_unique<InternalData>(QProjector<dim>::project_to_all_subfaces(
163 ReferenceCells::get_hypercube<dim>(), quadrature));
164 auto &data = dynamic_cast<InternalData &>(*data_ptr);
165
166 // verify that we have computed the transitive hull of the required
167 // flags and that FEValues has faithfully passed them on to us
168 Assert(update_flags == requires_update_flags(update_flags),
170
171 // store the flags in the internal data object so we can access them
172 // in fill_fe_*_values()
173 data.update_each = update_flags;
174
175 return data_ptr;
176}
177
178
179
180template <int dim, int spacedim>
181void
184 const CellSimilarity::Similarity cell_similarity,
185 const InternalData & data) const
186{
187 // Compute start point and sizes
188 // along axes. Strange vertex
189 // numbering makes this complicated
190 // again.
191 if (cell_similarity != CellSimilarity::translation)
192 {
193 const Point<dim> &start = cell->vertex(0);
194 switch (dim)
195 {
196 case 1:
197 data.cell_extents[0] = cell->vertex(1)(0) - start(0);
198 break;
199 case 2:
200 data.cell_extents[0] = cell->vertex(1)(0) - start(0);
201 data.cell_extents[1] = cell->vertex(2)(1) - start(1);
202 break;
203 case 3:
204 data.cell_extents[0] = cell->vertex(1)(0) - start(0);
205 data.cell_extents[1] = cell->vertex(2)(1) - start(1);
206 data.cell_extents[2] = cell->vertex(4)(2) - start(2);
207 break;
208 default:
209 Assert(false, ExcNotImplemented());
210 }
211 }
212}
213
214
215
216template <int dim, int spacedim>
217void
220 const InternalData & data,
221 std::vector<Point<dim>> &quadrature_points) const
222{
223 if (data.update_each & update_quadrature_points)
224 {
225 const auto offset = QProjector<dim>::DataSetDescriptor::cell();
226
228 }
229}
230
231
232
233template <int dim, int spacedim>
234void
237 const unsigned int face_no,
238 const InternalData & data,
239 std::vector<Point<dim>> &quadrature_points) const
240{
242
243 if (data.update_each & update_quadrature_points)
244 {
246 ReferenceCells::get_hypercube<dim>(),
247 face_no,
248 cell->face_orientation(face_no),
249 cell->face_flip(face_no),
250 cell->face_rotation(face_no),
251 quadrature_points.size());
252
253
255 }
256}
257
258
259
260template <int dim, int spacedim>
261void
264 const unsigned int face_no,
265 const unsigned int sub_no,
266 const InternalData & data,
267 std::vector<Point<dim>> &quadrature_points) const
268{
271 if (cell->face(face_no)->has_children())
272 {
273 AssertIndexRange(sub_no, cell->face(face_no)->n_children());
274 }
275
276 if (data.update_each & update_quadrature_points)
277 {
279 ReferenceCells::get_hypercube<dim>(),
280 face_no,
281 sub_no,
282 cell->face_orientation(face_no),
283 cell->face_flip(face_no),
284 cell->face_rotation(face_no),
285 quadrature_points.size(),
286 cell->subface_case(face_no));
287
289 }
290}
291
292
293
294template <int dim, int spacedim>
295void
298 const InternalData & data,
299 const typename QProjector<dim>::DataSetDescriptor & offset,
300 std::vector<Point<dim>> &quadrature_points) const
301{
302 // let @p{start} be the origin of a local coordinate system. it is chosen
303 // as the (lower) left vertex
304 const Point<dim> &start = cell->vertex(0);
305
306 for (unsigned int i = 0; i < quadrature_points.size(); ++i)
307 {
308 quadrature_points[i] = start;
309 for (unsigned int d = 0; d < dim; ++d)
310 quadrature_points[i](d) +=
311 data.cell_extents[d] * data.quadrature_points[i + offset](d);
312 }
313}
314
315
316
317template <int dim, int spacedim>
318void
320 const unsigned int face_no,
321 const InternalData & data,
322 std::vector<Tensor<1, dim>> &normal_vectors) const
323{
324 // compute normal vectors. All normals on a face have the same value.
325 if (data.update_each & update_normal_vectors)
326 {
328 std::fill(normal_vectors.begin(),
329 normal_vectors.end(),
331 }
332}
333
334
335
336template <int dim, int spacedim>
337void
339 const InternalData & data,
340 const CellSimilarity::Similarity cell_similarity,
342 &output_data) const
343{
344 if (cell_similarity != CellSimilarity::translation)
345 {
346 if (data.update_each & update_jacobian_grads)
347 for (unsigned int i = 0; i < output_data.jacobian_grads.size(); ++i)
349
350 if (data.update_each & update_jacobian_pushed_forward_grads)
351 for (unsigned int i = 0;
352 i < output_data.jacobian_pushed_forward_grads.size();
353 ++i)
355
356 if (data.update_each & update_jacobian_2nd_derivatives)
357 for (unsigned int i = 0;
358 i < output_data.jacobian_2nd_derivatives.size();
359 ++i)
360 output_data.jacobian_2nd_derivatives[i] =
362
364 for (unsigned int i = 0;
365 i < output_data.jacobian_pushed_forward_2nd_derivatives.size();
366 ++i)
369
370 if (data.update_each & update_jacobian_3rd_derivatives)
371 for (unsigned int i = 0;
372 i < output_data.jacobian_3rd_derivatives.size();
373 ++i)
374 output_data.jacobian_3rd_derivatives[i] =
376
378 for (unsigned int i = 0;
379 i < output_data.jacobian_pushed_forward_3rd_derivatives.size();
380 ++i)
383 }
384}
385
386
387
388template <int dim, int spacedim>
392 const CellSimilarity::Similarity cell_similarity,
393 const Quadrature<dim> & quadrature,
394 const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
396 &output_data) const
397{
398 // convert data object to internal data for this class. fails with
399 // an exception if that is not possible
400 Assert(dynamic_cast<const InternalData *>(&internal_data) != nullptr,
402 const InternalData &data = static_cast<const InternalData &>(internal_data);
403
404
405 update_cell_extents(cell, cell_similarity, data);
406
408 data,
409 output_data.quadrature_points);
410
411 // compute Jacobian determinant. all values are equal and are the
412 // product of the local lengths in each coordinate direction
413 if (data.update_each & (update_JxW_values | update_volume_elements))
414 if (cell_similarity != CellSimilarity::translation)
415 {
416 double J = data.cell_extents[0];
417 for (unsigned int d = 1; d < dim; ++d)
418 J *= data.cell_extents[d];
419 data.volume_element = J;
420 if (data.update_each & update_JxW_values)
421 for (unsigned int i = 0; i < output_data.JxW_values.size(); ++i)
422 output_data.JxW_values[i] = J * quadrature.weight(i);
423 }
424 // "compute" Jacobian at the quadrature points, which are all the
425 // same
426 if (data.update_each & update_jacobians)
427 if (cell_similarity != CellSimilarity::translation)
428 for (unsigned int i = 0; i < output_data.jacobians.size(); ++i)
429 {
431 for (unsigned int j = 0; j < dim; ++j)
432 output_data.jacobians[i][j][j] = data.cell_extents[j];
433 }
434
435 maybe_update_jacobian_derivatives(data, cell_similarity, output_data);
436
437 // "compute" inverse Jacobian at the quadrature points, which are
438 // all the same
439 if (data.update_each & update_inverse_jacobians)
440 if (cell_similarity != CellSimilarity::translation)
441 for (unsigned int i = 0; i < output_data.inverse_jacobians.size(); ++i)
442 {
443 output_data.inverse_jacobians[i] = Tensor<2, dim>();
444 for (unsigned int j = 0; j < dim; ++j)
445 output_data.inverse_jacobians[i][j][j] = 1. / data.cell_extents[j];
446 }
447
448 return cell_similarity;
449}
450
451
452
453template <int dim, int spacedim>
454void
457 const unsigned int face_no,
458 const hp::QCollection<dim - 1> & quadrature,
459 const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
461 &output_data) const
462{
463 AssertDimension(quadrature.size(), 1);
464
465 // convert data object to internal
466 // data for this class. fails with
467 // an exception if that is not
468 // possible
469 Assert(dynamic_cast<const InternalData *>(&internal_data) != nullptr,
471 const InternalData &data = static_cast<const InternalData &>(internal_data);
472
474
476 face_no,
477 data,
478 output_data.quadrature_points);
479
480 maybe_update_normal_vectors(face_no, data, output_data.normal_vectors);
481
482 // first compute Jacobian determinant, which is simply the product
483 // of the local lengths since the jacobian is diagonal
484 double J = 1.;
485 for (unsigned int d = 0; d < dim; ++d)
487 J *= data.cell_extents[d];
488
489 if (data.update_each & update_JxW_values)
490 for (unsigned int i = 0; i < output_data.JxW_values.size(); ++i)
491 output_data.JxW_values[i] = J * quadrature[0].weight(i);
492
493 if (data.update_each & update_boundary_forms)
494 for (unsigned int i = 0; i < output_data.boundary_forms.size(); ++i)
495 output_data.boundary_forms[i] = J * output_data.normal_vectors[i];
496
497 if (data.update_each & update_volume_elements)
498 {
499 J = data.cell_extents[0];
500 for (unsigned int d = 1; d < dim; ++d)
501 J *= data.cell_extents[d];
502 data.volume_element = J;
503 }
504
505 if (data.update_each & update_jacobians)
506 for (unsigned int i = 0; i < output_data.jacobians.size(); ++i)
507 {
509 for (unsigned int d = 0; d < dim; ++d)
510 output_data.jacobians[i][d][d] = data.cell_extents[d];
511 }
512
514
515 if (data.update_each & update_inverse_jacobians)
516 for (unsigned int i = 0; i < output_data.inverse_jacobians.size(); ++i)
517 {
519 for (unsigned int d = 0; d < dim; ++d)
520 output_data.inverse_jacobians[i][d][d] = 1. / data.cell_extents[d];
521 }
522}
523
524
525
526template <int dim, int spacedim>
527void
530 const unsigned int face_no,
531 const unsigned int subface_no,
532 const Quadrature<dim - 1> & quadrature,
533 const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
535 &output_data) const
536{
537 // convert data object to internal data for this class. fails with
538 // an exception if that is not possible
539 Assert(dynamic_cast<const InternalData *>(&internal_data) != nullptr,
541 const InternalData &data = static_cast<const InternalData &>(internal_data);
542
544
546 cell, face_no, subface_no, data, output_data.quadrature_points);
547
548 maybe_update_normal_vectors(face_no, data, output_data.normal_vectors);
549
550 // first compute Jacobian determinant, which is simply the product
551 // of the local lengths since the jacobian is diagonal
552 double J = 1.;
553 for (unsigned int d = 0; d < dim; ++d)
555 J *= data.cell_extents[d];
556
557 if (data.update_each & update_JxW_values)
558 {
559 // Here, cell->face(face_no)->n_children() would be the right
560 // choice, but unfortunately the current function is also called
561 // for faces without children (see tests/fe/mapping.cc). Add
562 // following switch to avoid diffs in tests/fe/mapping.OK
563 const unsigned int n_subfaces =
564 cell->face(face_no)->has_children() ?
565 cell->face(face_no)->n_children() :
567 for (unsigned int i = 0; i < output_data.JxW_values.size(); ++i)
568 output_data.JxW_values[i] = J * quadrature.weight(i) / n_subfaces;
569 }
570
571 if (data.update_each & update_boundary_forms)
572 for (unsigned int i = 0; i < output_data.boundary_forms.size(); ++i)
573 output_data.boundary_forms[i] = J * output_data.normal_vectors[i];
574
575 if (data.update_each & update_volume_elements)
576 {
577 J = data.cell_extents[0];
578 for (unsigned int d = 1; d < dim; ++d)
579 J *= data.cell_extents[d];
580 data.volume_element = J;
581 }
582
583 if (data.update_each & update_jacobians)
584 for (unsigned int i = 0; i < output_data.jacobians.size(); ++i)
585 {
587 for (unsigned int d = 0; d < dim; ++d)
588 output_data.jacobians[i][d][d] = data.cell_extents[d];
589 }
590
592
593 if (data.update_each & update_inverse_jacobians)
594 for (unsigned int i = 0; i < output_data.inverse_jacobians.size(); ++i)
595 {
597 for (unsigned int d = 0; d < dim; ++d)
598 output_data.inverse_jacobians[i][d][d] = 1. / data.cell_extents[d];
599 }
600}
601
602
603
604template <int dim, int spacedim>
605void
607 const ArrayView<const Tensor<1, dim>> & input,
608 const MappingKind mapping_kind,
609 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
610 const ArrayView<Tensor<1, spacedim>> & output) const
611{
612 AssertDimension(input.size(), output.size());
613 Assert(dynamic_cast<const InternalData *>(&mapping_data) != nullptr,
615 const InternalData &data = static_cast<const InternalData &>(mapping_data);
616
617 switch (mapping_kind)
618 {
620 {
621 Assert(data.update_each & update_covariant_transformation,
623 "update_covariant_transformation"));
624
625 for (unsigned int i = 0; i < output.size(); ++i)
626 for (unsigned int d = 0; d < dim; ++d)
627 output[i][d] = input[i][d] / data.cell_extents[d];
628 return;
629 }
630
632 {
635 "update_contravariant_transformation"));
636
637 for (unsigned int i = 0; i < output.size(); ++i)
638 for (unsigned int d = 0; d < dim; ++d)
639 output[i][d] = input[i][d] * data.cell_extents[d];
640 return;
641 }
642 case mapping_piola:
643 {
646 "update_contravariant_transformation"));
647 Assert(data.update_each & update_volume_elements,
649 "update_volume_elements"));
650
651 for (unsigned int i = 0; i < output.size(); ++i)
652 for (unsigned int d = 0; d < dim; ++d)
653 output[i][d] =
654 input[i][d] * data.cell_extents[d] / data.volume_element;
655 return;
656 }
657 default:
658 Assert(false, ExcNotImplemented());
659 }
660}
661
662
663
664template <int dim, int spacedim>
665void
668 const MappingKind mapping_kind,
669 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
670 const ArrayView<Tensor<2, spacedim>> & output) const
671{
672 AssertDimension(input.size(), output.size());
673 Assert(dynamic_cast<const InternalData *>(&mapping_data) != nullptr,
675 const InternalData &data = static_cast<const InternalData &>(mapping_data);
676
677 switch (mapping_kind)
678 {
680 {
681 Assert(data.update_each & update_covariant_transformation,
683 "update_covariant_transformation"));
684
685 for (unsigned int i = 0; i < output.size(); ++i)
686 for (unsigned int d1 = 0; d1 < dim; ++d1)
687 for (unsigned int d2 = 0; d2 < dim; ++d2)
688 output[i][d1][d2] = input[i][d1][d2] / data.cell_extents[d2];
689 return;
690 }
691
693 {
696 "update_contravariant_transformation"));
697
698 for (unsigned int i = 0; i < output.size(); ++i)
699 for (unsigned int d1 = 0; d1 < dim; ++d1)
700 for (unsigned int d2 = 0; d2 < dim; ++d2)
701 output[i][d1][d2] = input[i][d1][d2] * data.cell_extents[d2];
702 return;
703 }
704
706 {
707 Assert(data.update_each & update_covariant_transformation,
709 "update_covariant_transformation"));
710
711 for (unsigned int i = 0; i < output.size(); ++i)
712 for (unsigned int d1 = 0; d1 < dim; ++d1)
713 for (unsigned int d2 = 0; d2 < dim; ++d2)
714 output[i][d1][d2] = input[i][d1][d2] / data.cell_extents[d2] /
715 data.cell_extents[d1];
716 return;
717 }
718
720 {
723 "update_contravariant_transformation"));
724
725 for (unsigned int i = 0; i < output.size(); ++i)
726 for (unsigned int d1 = 0; d1 < dim; ++d1)
727 for (unsigned int d2 = 0; d2 < dim; ++d2)
728 output[i][d1][d2] = input[i][d1][d2] * data.cell_extents[d2] /
729 data.cell_extents[d1];
730 return;
731 }
732
733 case mapping_piola:
734 {
737 "update_contravariant_transformation"));
738 Assert(data.update_each & update_volume_elements,
740 "update_volume_elements"));
741
742 for (unsigned int i = 0; i < output.size(); ++i)
743 for (unsigned int d1 = 0; d1 < dim; ++d1)
744 for (unsigned int d2 = 0; d2 < dim; ++d2)
745 output[i][d1][d2] = input[i][d1][d2] * data.cell_extents[d2] /
746 data.volume_element;
747 return;
748 }
749
751 {
754 "update_contravariant_transformation"));
755 Assert(data.update_each & update_volume_elements,
757 "update_volume_elements"));
758
759 for (unsigned int i = 0; i < output.size(); ++i)
760 for (unsigned int d1 = 0; d1 < dim; ++d1)
761 for (unsigned int d2 = 0; d2 < dim; ++d2)
762 output[i][d1][d2] = input[i][d1][d2] * data.cell_extents[d2] /
763 data.cell_extents[d1] / data.volume_element;
764 return;
765 }
766
767 default:
768 Assert(false, ExcNotImplemented());
769 }
770}
771
772
773
774template <int dim, int spacedim>
775void
777 const ArrayView<const Tensor<2, dim>> & input,
778 const MappingKind mapping_kind,
779 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
780 const ArrayView<Tensor<2, spacedim>> & output) const
781{
782 AssertDimension(input.size(), output.size());
783 Assert(dynamic_cast<const InternalData *>(&mapping_data) != nullptr,
785 const InternalData &data = static_cast<const InternalData &>(mapping_data);
786
787 switch (mapping_kind)
788 {
790 {
791 Assert(data.update_each & update_covariant_transformation,
793 "update_covariant_transformation"));
794
795 for (unsigned int i = 0; i < output.size(); ++i)
796 for (unsigned int d1 = 0; d1 < dim; ++d1)
797 for (unsigned int d2 = 0; d2 < dim; ++d2)
798 output[i][d1][d2] = input[i][d1][d2] / data.cell_extents[d2];
799 return;
800 }
801
803 {
806 "update_contravariant_transformation"));
807
808 for (unsigned int i = 0; i < output.size(); ++i)
809 for (unsigned int d1 = 0; d1 < dim; ++d1)
810 for (unsigned int d2 = 0; d2 < dim; ++d2)
811 output[i][d1][d2] = input[i][d1][d2] * data.cell_extents[d2];
812 return;
813 }
814
816 {
817 Assert(data.update_each & update_covariant_transformation,
819 "update_covariant_transformation"));
820
821 for (unsigned int i = 0; i < output.size(); ++i)
822 for (unsigned int d1 = 0; d1 < dim; ++d1)
823 for (unsigned int d2 = 0; d2 < dim; ++d2)
824 output[i][d1][d2] = input[i][d1][d2] / data.cell_extents[d2] /
825 data.cell_extents[d1];
826 return;
827 }
828
830 {
833 "update_contravariant_transformation"));
834
835 for (unsigned int i = 0; i < output.size(); ++i)
836 for (unsigned int d1 = 0; d1 < dim; ++d1)
837 for (unsigned int d2 = 0; d2 < dim; ++d2)
838 output[i][d1][d2] = input[i][d1][d2] * data.cell_extents[d2] /
839 data.cell_extents[d1];
840 return;
841 }
842
843 case mapping_piola:
844 {
847 "update_contravariant_transformation"));
848 Assert(data.update_each & update_volume_elements,
850 "update_volume_elements"));
851
852 for (unsigned int i = 0; i < output.size(); ++i)
853 for (unsigned int d1 = 0; d1 < dim; ++d1)
854 for (unsigned int d2 = 0; d2 < dim; ++d2)
855 output[i][d1][d2] = input[i][d1][d2] * data.cell_extents[d2] /
856 data.volume_element;
857 return;
858 }
859
861 {
864 "update_contravariant_transformation"));
865 Assert(data.update_each & update_volume_elements,
867 "update_volume_elements"));
868
869 for (unsigned int i = 0; i < output.size(); ++i)
870 for (unsigned int d1 = 0; d1 < dim; ++d1)
871 for (unsigned int d2 = 0; d2 < dim; ++d2)
872 output[i][d1][d2] = input[i][d1][d2] * data.cell_extents[d2] /
873 data.cell_extents[d1] / data.volume_element;
874 return;
875 }
876
877 default:
878 Assert(false, ExcNotImplemented());
879 }
880}
881
882
883
884template <int dim, int spacedim>
885void
888 const MappingKind mapping_kind,
889 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
890 const ArrayView<Tensor<3, spacedim>> & output) const
891{
892 AssertDimension(input.size(), output.size());
893 Assert(dynamic_cast<const InternalData *>(&mapping_data) != nullptr,
895 const InternalData &data = static_cast<const InternalData &>(mapping_data);
896
897 switch (mapping_kind)
898 {
900 {
903 "update_covariant_transformation"));
904
905 for (unsigned int q = 0; q < output.size(); ++q)
906 for (unsigned int i = 0; i < spacedim; ++i)
907 for (unsigned int j = 0; j < spacedim; ++j)
908 for (unsigned int k = 0; k < spacedim; ++k)
909 {
910 output[q][i][j][k] = input[q][i][j][k] /
911 data.cell_extents[j] /
912 data.cell_extents[k];
913 }
914 return;
915 }
916 default:
917 Assert(false, ExcNotImplemented());
918 }
919}
920
921
922
923template <int dim, int spacedim>
924void
926 const ArrayView<const Tensor<3, dim>> & input,
927 const MappingKind mapping_kind,
928 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
929 const ArrayView<Tensor<3, spacedim>> & output) const
930{
931 AssertDimension(input.size(), output.size());
932 Assert(dynamic_cast<const InternalData *>(&mapping_data) != nullptr,
934 const InternalData &data = static_cast<const InternalData &>(mapping_data);
935
936 switch (mapping_kind)
937 {
939 {
940 Assert(data.update_each & update_covariant_transformation,
942 "update_covariant_transformation"));
945 "update_contravariant_transformation"));
946
947 for (unsigned int q = 0; q < output.size(); ++q)
948 for (unsigned int i = 0; i < spacedim; ++i)
949 for (unsigned int j = 0; j < spacedim; ++j)
950 for (unsigned int k = 0; k < spacedim; ++k)
951 {
952 output[q][i][j][k] =
953 input[q][i][j][k] * data.cell_extents[i] /
954 data.cell_extents[j] / data.cell_extents[k];
955 }
956 return;
957 }
958
960 {
961 Assert(data.update_each & update_covariant_transformation,
963 "update_covariant_transformation"));
964
965 for (unsigned int q = 0; q < output.size(); ++q)
966 for (unsigned int i = 0; i < spacedim; ++i)
967 for (unsigned int j = 0; j < spacedim; ++j)
968 for (unsigned int k = 0; k < spacedim; ++k)
969 {
970 output[q][i][j][k] =
971 input[q][i][j][k] / data.cell_extents[i] /
972 data.cell_extents[j] / data.cell_extents[k];
973 }
974
975 return;
976 }
977
979 {
980 Assert(data.update_each & update_covariant_transformation,
982 "update_covariant_transformation"));
985 "update_contravariant_transformation"));
986 Assert(data.update_each & update_volume_elements,
988 "update_volume_elements"));
989
990 for (unsigned int q = 0; q < output.size(); ++q)
991 for (unsigned int i = 0; i < spacedim; ++i)
992 for (unsigned int j = 0; j < spacedim; ++j)
993 for (unsigned int k = 0; k < spacedim; ++k)
994 {
995 output[q][i][j][k] =
996 input[q][i][j][k] * data.cell_extents[i] /
997 data.volume_element / data.cell_extents[j] /
998 data.cell_extents[k];
999 }
1000
1001 return;
1002 }
1003
1004 default:
1005 Assert(false, ExcNotImplemented());
1006 }
1007}
1008
1009
1010
1011template <int dim, int spacedim>
1015 const Point<dim> & p) const
1016{
1017 Tensor<1, dim> length;
1018 const Point<dim> start = cell->vertex(0);
1019 switch (dim)
1020 {
1021 case 1:
1022 length[0] = cell->vertex(1)(0) - start(0);
1023 break;
1024 case 2:
1025 length[0] = cell->vertex(1)(0) - start(0);
1026 length[1] = cell->vertex(2)(1) - start(1);
1027 break;
1028 case 3:
1029 length[0] = cell->vertex(1)(0) - start(0);
1030 length[1] = cell->vertex(2)(1) - start(1);
1031 length[2] = cell->vertex(4)(2) - start(2);
1032 break;
1033 default:
1034 Assert(false, ExcNotImplemented());
1035 }
1036
1037 Point<dim> p_real = cell->vertex(0);
1038 for (unsigned int d = 0; d < dim; ++d)
1039 p_real(d) += length[d] * p(d);
1040
1041 return p_real;
1042}
1043
1044
1045
1046template <int dim, int spacedim>
1050 const Point<spacedim> & p) const
1051{
1052 if (dim != spacedim)
1053 Assert(false, ExcNotImplemented());
1054 const Point<dim> &start = cell->vertex(0);
1055 Point<dim> real = p;
1056 real -= start;
1057
1058 switch (dim)
1059 {
1060 case 1:
1061 real(0) /= cell->vertex(1)(0) - start(0);
1062 break;
1063 case 2:
1064 real(0) /= cell->vertex(1)(0) - start(0);
1065 real(1) /= cell->vertex(2)(1) - start(1);
1066 break;
1067 case 3:
1068 real(0) /= cell->vertex(1)(0) - start(0);
1069 real(1) /= cell->vertex(2)(1) - start(1);
1070 real(2) /= cell->vertex(4)(2) - start(2);
1071 break;
1072 default:
1073 Assert(false, ExcNotImplemented());
1074 }
1075 return real;
1076}
1077
1078
1079
1080template <int dim, int spacedim>
1081std::unique_ptr<Mapping<dim, spacedim>>
1083{
1084 return std::make_unique<MappingCartesian<dim, spacedim>>(*this);
1085}
1086
1087
1088//---------------------------------------------------------------------------
1089// explicit instantiations
1090#include "mapping_cartesian.inst"
1091
1092
std::vector< Point< dim > > quadrature_points
InternalData(const Quadrature< dim > &quadrature)
virtual std::size_t memory_consumption() const override
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const hp::QCollection< dim - 1 > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_face_data(const UpdateFlags flags, const hp::QCollection< dim - 1 > &quadrature) const override
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags flags, const Quadrature< dim - 1 > &quadrature) const override
void update_cell_extents(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const InternalData &data) const
void maybe_update_subface_quadrature_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const InternalData &data, std::vector< Point< dim > > &quadrature_points) const
virtual CellSimilarity::Similarity fill_fe_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const Quadrature< dim > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
void maybe_update_face_quadrature_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const InternalData &data, std::vector< Point< dim > > &quadrature_points) const
virtual bool is_compatible_with(const ReferenceCell &reference_cell) const override
void maybe_update_cell_quadrature_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const InternalData &data, std::vector< Point< dim > > &quadrature_points) const
virtual bool preserves_vertex_locations() const override
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const override
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const override
void transform_quadrature_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const InternalData &data, const typename QProjector< dim >::DataSetDescriptor &offset, std::vector< Point< dim > > &quadrature_points) const
void maybe_update_jacobian_derivatives(const InternalData &data, const CellSimilarity::Similarity cell_similarity, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const
virtual void transform(const ArrayView< const Tensor< 1, dim > > &input, const MappingKind kind, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const ArrayView< Tensor< 1, spacedim > > &output) const override
void maybe_update_normal_vectors(const unsigned int face_no, const InternalData &data, std::vector< Tensor< 1, dim > > &normal_vectors) const
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const override
virtual void fill_fe_subface_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int subface_no, const Quadrature< dim - 1 > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
Abstract base class for mapping classes.
Definition: mapping.h:304
Definition: point.h:111
static DataSetDescriptor cell()
Definition: qprojector.h:563
static DataSetDescriptor subface(const unsigned int face_no, const unsigned int subface_no, const bool face_orientation, const bool face_flip, const bool face_rotation, const unsigned int n_quadrature_points, const internal::SubfaceCase< dim > ref_case=internal::SubfaceCase< dim >::case_isotropic)
static DataSetDescriptor face(const unsigned int face_no, const bool face_orientation, const bool face_flip, const bool face_rotation, const unsigned int n_quadrature_points)
Definition: qprojector.cc:1365
double weight(const unsigned int i) const
Definition: tensor.h:472
unsigned int size() const
Definition: collection.h:109
std::vector< DerivativeForm< 1, spacedim, dim > > inverse_jacobians
std::vector< Tensor< 1, spacedim > > normal_vectors
std::vector< Tensor< 5, spacedim > > jacobian_pushed_forward_3rd_derivatives
std::vector< DerivativeForm< 4, dim, spacedim > > jacobian_3rd_derivatives
std::vector< DerivativeForm< 3, dim, spacedim > > jacobian_2nd_derivatives
std::vector< Point< spacedim > > quadrature_points
std::vector< Tensor< 4, spacedim > > jacobian_pushed_forward_2nd_derivatives
std::vector< Tensor< 3, spacedim > > jacobian_pushed_forward_grads
std::vector< DerivativeForm< 2, dim, spacedim > > jacobian_grads
std::vector< Tensor< 1, spacedim > > boundary_forms
std::vector< DerivativeForm< 1, dim, spacedim > > jacobians
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:402
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:403
UpdateFlags
@ update_jacobian_pushed_forward_2nd_derivatives
@ update_volume_elements
Determinant of the Jacobian.
@ update_contravariant_transformation
Contravariant transformation.
@ update_jacobian_pushed_forward_grads
@ update_jacobian_3rd_derivatives
@ update_jacobian_grads
Gradient of volume element.
@ update_normal_vectors
Normal vectors.
@ update_JxW_values
Transformed quadrature weights.
@ update_covariant_transformation
Covariant transformation.
@ update_jacobians
Volume element.
@ update_inverse_jacobians
Volume element.
@ update_quadrature_points
Transformed quadrature points.
@ update_jacobian_pushed_forward_3rd_derivatives
@ update_boundary_forms
Outer normal vector, not normalized.
@ update_jacobian_2nd_derivatives
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
Definition: exceptions.h:1465
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1622
#define AssertIndexRange(index, range)
Definition: exceptions.h:1690
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcMessage(std::string arg1)
MappingKind
Definition: mapping.h:65
@ mapping_piola
Definition: mapping.h:100
@ mapping_covariant_gradient
Definition: mapping.h:86
@ mapping_covariant
Definition: mapping.h:75
@ mapping_contravariant
Definition: mapping.h:80
@ mapping_contravariant_hessian
Definition: mapping.h:142
@ mapping_covariant_hessian
Definition: mapping.h:136
@ mapping_contravariant_gradient
Definition: mapping.h:92
@ mapping_piola_gradient
Definition: mapping.h:106
@ mapping_piola_hessian
Definition: mapping.h:148
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
void quadrature_points(const Triangulation< dim, spacedim > &triangulation, const Quadrature< dim > &quadrature, const std::vector< std::vector< BoundingBox< spacedim > > > &global_bounding_boxes, ParticleHandler< dim, spacedim > &particle_handler, const Mapping< dim, spacedim > &mapping=(ReferenceCells::get_hypercube< dim >() .template get_default_linear_mapping< dim, spacedim >()), const std::vector< std::vector< double > > &properties={})
Definition: generators.cc:451
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
std::string to_string(const number value, const unsigned int digits=numbers::invalid_unsigned_int)
Definition: utilities.cc:482
T signaling_nan()