Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
mapping_info.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2011 - 2023 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16
17#ifndef dealii_non_matching_mapping_info_h
18#define dealii_non_matching_mapping_info_h
19
20
21#include <deal.II/base/config.h>
22
26
27#include <deal.II/fe/fe_dgq.h>
29#include <deal.II/fe/mapping.h>
33
34#include <memory>
35
36
38
39namespace NonMatching
40{
41 namespace internal
42 {
43 template <int dim, int spacedim = dim>
45 {
48 spacedim>;
49
50 public:
51 static UpdateFlags
53 const SmartPointer<const Mapping<dim, spacedim>> mapping,
54 const UpdateFlags & update_flags)
55 {
56 return mapping->requires_update_flags(update_flags);
57 }
58
59 static void
61 const SmartPointer<const Mapping<dim, spacedim>> mapping,
62 const UpdateFlags & update_flags_mapping,
64 CellSimilarity::Similarity &cell_similarity,
65 const Quadrature<dim> & quadrature,
66 std::shared_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
67 internal_mapping_data,
68 MappingData &mapping_data)
69 {
70 mapping_data.initialize(quadrature.size(), update_flags_mapping);
71
72 // reuse internal_mapping_data for MappingQ to avoid memory allocations
73 if (const MappingQ<dim, spacedim> *mapping_q =
74 dynamic_cast<const MappingQ<dim, spacedim> *>(&(*mapping)))
75 {
76 (void)mapping_q;
77 auto &data =
78 dynamic_cast<typename MappingQ<dim, spacedim>::InternalData &>(
79 *internal_mapping_data);
80 data.initialize(update_flags_mapping,
81 quadrature,
82 quadrature.size());
83 }
84 else
85 {
86 internal_mapping_data =
87 mapping->get_data(update_flags_mapping, quadrature);
88 }
89
90 cell_similarity = mapping->fill_fe_values(cell,
91 cell_similarity,
92 quadrature,
93 *internal_mapping_data,
94 mapping_data);
95 }
96
97
98
99 static void
101 const SmartPointer<const Mapping<dim, spacedim>> mapping,
102 const UpdateFlags & update_flags_mapping,
104 const ImmersedSurfaceQuadrature<dim> & quadrature,
105 std::shared_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
106 internal_mapping_data,
107 MappingData &mapping_data)
108 {
109 mapping_data.initialize(quadrature.size(), update_flags_mapping);
110
111 // reuse internal_mapping_data for MappingQ to avoid memory allocations
112 if (dynamic_cast<const MappingQ<dim, spacedim> *>(&(*mapping)))
113 {
114 auto &data =
115 dynamic_cast<typename MappingQ<dim, spacedim>::InternalData &>(
116 *internal_mapping_data);
117 data.initialize(update_flags_mapping,
118 quadrature,
119 quadrature.size());
120 }
121 else
122 {
123 internal_mapping_data =
124 mapping->get_data(update_flags_mapping, quadrature);
125 }
126
127 mapping->fill_fe_immersed_surface_values(cell,
128 quadrature,
129 *internal_mapping_data,
130 mapping_data);
131 }
132
133
134
135 static void
137 const SmartPointer<const Mapping<dim, spacedim>> mapping,
138 const UpdateFlags & update_flags_mapping,
140 const unsigned int face_no,
141 const Quadrature<dim - 1> & quadrature,
142 std::shared_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
143 internal_mapping_data,
144 MappingData &mapping_data)
145 {
146 mapping_data.initialize(quadrature.size(), update_flags_mapping);
147
148 // reuse internal_mapping_data for MappingQ to avoid memory allocations
149 if (const MappingQ<dim, spacedim> *mapping_q =
150 dynamic_cast<const MappingQ<dim, spacedim> *>(&(*mapping)))
151 {
152 auto &data =
153 dynamic_cast<typename MappingQ<dim, spacedim>::InternalData &>(
154 *internal_mapping_data);
155 data.initialize_face(update_flags_mapping,
157 ReferenceCells::get_hypercube<dim>(),
158 quadrature,
159 face_no,
160 cell->face_orientation(face_no),
161 cell->face_flip(face_no),
162 cell->face_rotation(face_no)),
163 quadrature.size());
164
165 mapping_q->fill_mapping_data_for_face_quadrature(
166 cell, face_no, quadrature, *internal_mapping_data, mapping_data);
167 }
168 else
169 {
170 auto internal_mapping_data =
171 mapping->get_face_data(update_flags_mapping,
172 hp::QCollection<dim - 1>(quadrature));
173
174 mapping->fill_fe_face_values(cell,
175 face_no,
176 hp::QCollection<dim - 1>(quadrature),
177 *internal_mapping_data,
178 mapping_data);
179 }
180 }
181 };
182 } // namespace internal
183
200 template <int dim, int spacedim = dim, typename Number = double>
202 {
203 public:
207 using VectorizedArrayType = typename ::internal::VectorizedArrayTrait<
208 Number>::vectorized_value_type;
209
215 {
221 {}
222
229 };
230
245 const AdditionalData additional_data = AdditionalData());
246
251 void
253 const std::vector<Point<dim>> &unit_points);
254
259 void
261 const ArrayView<const Point<dim>> &unit_points);
262
268 void
270 const Quadrature<dim> &quadrature);
271
282 template <typename ContainerType>
283 void
285 const ContainerType & cell_iterator_range,
286 const std::vector<std::vector<Point<dim>>> &unit_points_vector,
287 const unsigned int n_unfiltered_cells = numbers::invalid_unsigned_int);
288
295 template <typename ContainerType>
296 void
298 const ContainerType & cell_iterator_range,
299 const std::vector<Quadrature<dim>> &quadrature_vector,
300 const unsigned int n_unfiltered_cells = numbers::invalid_unsigned_int);
301
306 template <typename Iterator>
307 void
309 const IteratorRange<Iterator> & cell_iterator_range,
310 const std::vector<ImmersedSurfaceQuadrature<dim>> &quadrature_vector,
311 const unsigned int n_unfiltered_cells = numbers::invalid_unsigned_int);
312
317 template <typename Iterator>
318 void
320 const IteratorRange<Iterator> & cell_iterator_range,
321 const std::vector<std::vector<Quadrature<dim - 1>>> &quadrature_vector,
322 const unsigned int n_unfiltered_cells = numbers::invalid_unsigned_int);
323
328 bool
330
336 get_unit_point(const unsigned int offset) const;
337
342 const Point<dim - 1, VectorizedArrayType> *
343 get_unit_point_faces(const unsigned int offset) const;
344
350 get_jacobian(const unsigned int offset) const;
351
357 get_inverse_jacobian(const unsigned int offset) const;
358
364 get_normal_vector(const unsigned int offset) const;
365
370 const Number *
371 get_JxW(const unsigned int offset) const;
372
377 const Point<dim, Number> *
378 get_real_point(const unsigned int offset) const;
379
384 get_mapping() const;
385
391
397
401 boost::signals2::connection
402 connect_is_reinitialized(const std::function<void()> &set_is_reinitialized);
403
407 unsigned int
409 const unsigned int face_number) const;
410
414 unsigned int
416 const unsigned int face_number) const;
417
421 unsigned int
423 const unsigned int face_number) const;
424
425 private:
428 spacedim>;
429
433 template <typename NumberType>
434 unsigned int
436
440 void
442
446 void
447 resize_unit_points(const unsigned int n_unit_point_batches);
448
452 void
453 resize_unit_points_faces(const unsigned int n_unit_point_batches);
454
458 void
459 resize_data_fields(const unsigned int n_data_point_batches);
460
464 void
465 store_unit_points(const unsigned int unit_points_index_offset,
466 const unsigned int n_q_points,
467 const unsigned int n_q_points_unvectorized,
468 const std::vector<Point<dim>> &points);
469
473 void
474 store_unit_points_faces(const unsigned int unit_points_index_offset,
475 const unsigned int n_q_points,
476 const unsigned int n_q_points_unvectorized,
477 const std::vector<Point<dim - 1>> &points);
478
482 void
483 store_mapping_data(const unsigned int unit_points_index_offset,
484 const unsigned int n_q_points,
485 const unsigned int n_q_points_unvectorized,
486 const MappingData & mapping_data,
487 const std::vector<double> &weights);
488
492 unsigned int
493 compute_compressed_cell_index(const unsigned int cell_index) const;
494
498 unsigned int
500 const unsigned int face_number) const;
501
505 enum class State
506 {
507 invalid,
511 };
512
518
525
532
537
541 std::shared_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
543
548
553
558
563
569
577
584
592
601
608
612 std::vector<unsigned int> n_q_points_unvectorized;
613
618 std::vector<unsigned int> cell_index_offset;
619
624 std::vector<unsigned int> cell_index_to_compressed_cell_index;
625
630
635 boost::signals2::signal<void()> is_reinitialized;
636 };
637
638 // ----------------------- template functions ----------------------
639
640
641 template <int dim, int spacedim, typename Number>
643 const Mapping<dim> & mapping,
644 const UpdateFlags update_flags,
645 const AdditionalData additional_data)
646 : mapping(&mapping)
647 , update_flags(update_flags)
648 , update_flags_mapping(update_default)
649 , additional_data(additional_data)
650 {
651 // translate update flags
661
662 // always save quadrature points for now
664
667 this->mapping, update_flags_mapping);
668
669 // construct internal_mapping_data for MappingQ to be able to reuse it in
670 // reinit() calls to avoid memory allocations
671 if (const MappingQ<dim, spacedim> *mapping_q =
672 dynamic_cast<const MappingQ<dim, spacedim> *>(&mapping))
673 {
675 std::make_unique<typename MappingQ<dim, spacedim>::InternalData>(
676 mapping_q->get_degree());
677 }
678 }
679
680
681
682 template <int dim, int spacedim, typename Number>
683 void
685 {
686 n_q_points_unvectorized.clear();
687 unit_points_index.clear();
688 data_index_offsets.clear();
689 }
690
691
692
693 template <int dim, int spacedim, typename Number>
694 void
697 const std::vector<Point<dim>> & unit_points_in)
698 {
699 reinit(cell, Quadrature<dim>(unit_points_in));
700 }
701
702
703
704 template <int dim, int spacedim, typename Number>
705 void
708 const ArrayView<const Point<dim>> & unit_points_in)
709 {
710 reinit(cell,
711 std::vector<Point<dim>>(unit_points_in.begin(),
712 unit_points_in.end()));
713 }
714
715
716
717 template <int dim, int spacedim, typename Number>
718 void
721 const Quadrature<dim> & quadrature)
722 {
723 n_q_points_unvectorized.resize(1);
724 n_q_points_unvectorized[0] = quadrature.size();
725
726 const unsigned int n_q_points =
727 compute_n_q_points<VectorizedArrayType>(n_q_points_unvectorized[0]);
728
729 const unsigned int n_q_points_data =
730 compute_n_q_points<Number>(n_q_points_unvectorized[0]);
731
732 // resize data vectors
733 resize_unit_points(n_q_points);
734 resize_data_fields(n_q_points_data);
735
736 // store unit points
737 store_unit_points(0,
738 n_q_points,
739 n_q_points_unvectorized[0],
740 quadrature.get_points());
741
742 // compute mapping data
743 MappingData mapping_data;
747 update_flags_mapping,
748 cell,
749 cell_similarity,
750 quadrature,
751 internal_mapping_data,
752 mapping_data);
753
754 // store mapping data
755 store_mapping_data(0,
756 n_q_points_data,
757 n_q_points_unvectorized[0],
758 mapping_data,
759 quadrature.get_weights());
760
761 state = State::single_cell;
762 is_reinitialized();
763 }
764
765
766
767 template <int dim, int spacedim, typename Number>
768 template <typename ContainerType>
769 void
771 const ContainerType & cell_iterator_range,
772 const std::vector<std::vector<Point<dim>>> &unit_points_vector,
773 const unsigned int n_unfiltered_cells)
774 {
775 const unsigned int n_cells = unit_points_vector.size();
776 AssertDimension(n_cells,
777 std::distance(cell_iterator_range.begin(),
778 cell_iterator_range.end()));
779
780 std::vector<Quadrature<dim>> quadrature_vector(n_cells);
781 for (unsigned int cell_index = 0; cell_index < n_cells; ++cell_index)
782 quadrature_vector[cell_index] =
783 Quadrature<dim>(unit_points_vector[cell_index]);
784
785 reinit_cells(cell_iterator_range, quadrature_vector, n_unfiltered_cells);
786 }
787
788
789
790 template <int dim, int spacedim, typename Number>
791 template <typename ContainerType>
792 void
794 const ContainerType & cell_iterator_range,
795 const std::vector<Quadrature<dim>> &quadrature_vector,
796 const unsigned int n_unfiltered_cells)
797 {
798 clear();
799
800 do_cell_index_compression =
801 n_unfiltered_cells != numbers::invalid_unsigned_int;
802
803 const unsigned int n_cells = quadrature_vector.size();
804 AssertDimension(n_cells,
805 std::distance(cell_iterator_range.begin(),
806 cell_iterator_range.end()));
807
808 n_q_points_unvectorized.reserve(n_cells);
809
810 // fill unit points index offset vector
811 unit_points_index.reserve(n_cells + 1);
812 unit_points_index.push_back(0);
813 data_index_offsets.reserve(n_cells + 1);
814 data_index_offsets.push_back(0);
815 for (const auto &quadrature : quadrature_vector)
816 {
817 const unsigned int n_points = quadrature.size();
818 n_q_points_unvectorized.push_back(n_points);
819
820 const unsigned int n_q_points =
821 compute_n_q_points<VectorizedArrayType>(n_points);
822 unit_points_index.push_back(unit_points_index.back() + n_q_points);
823
824 const unsigned int n_q_points_data =
825 compute_n_q_points<Number>(n_points);
826 data_index_offsets.push_back(data_index_offsets.back() +
827 n_q_points_data);
828 }
829
830 const unsigned int n_unit_points = unit_points_index.back();
831 const unsigned int n_data_points = data_index_offsets.back();
832
833 // resize data vectors
834 resize_unit_points(n_unit_points);
835 resize_data_fields(n_data_points);
836
837 if (do_cell_index_compression)
838 cell_index_to_compressed_cell_index.resize(n_unfiltered_cells,
840
841 MappingData mapping_data;
842 CellSimilarity::Similarity cell_similarity =
844 unsigned int cell_index = 0;
845 for (const auto &cell : cell_iterator_range)
846 {
847 // store unit points
848 const unsigned int n_q_points = compute_n_q_points<VectorizedArrayType>(
849 n_q_points_unvectorized[cell_index]);
850 store_unit_points(unit_points_index[cell_index],
851 n_q_points,
852 n_q_points_unvectorized[cell_index],
853 quadrature_vector[cell_index].get_points());
854
855 // compute mapping data
858 update_flags_mapping,
859 cell,
860 cell_similarity,
861 quadrature_vector[cell_index],
862 internal_mapping_data,
863 mapping_data);
864
865 // store mapping data
866 const unsigned int n_q_points_data =
867 compute_n_q_points<Number>(n_q_points_unvectorized[cell_index]);
868 store_mapping_data(data_index_offsets[cell_index],
869 n_q_points_data,
870 n_q_points_unvectorized[cell_index],
871 mapping_data,
872 quadrature_vector[cell_index].get_weights());
873
874 if (do_cell_index_compression)
875 cell_index_to_compressed_cell_index[cell->active_cell_index()] =
877
878 ++cell_index;
879 }
880
881 state = State::cell_vector;
882 is_reinitialized();
883 }
884
885
886
887 template <int dim, int spacedim, typename Number>
888 template <typename Iterator>
889 void
891 const IteratorRange<Iterator> & cell_iterator_range,
892 const std::vector<ImmersedSurfaceQuadrature<dim>> &quadrature_vector,
893 const unsigned int n_unfiltered_cells)
894 {
895 clear();
896
897 Assert(
898 additional_data.use_global_weights == false,
900 "There is no known use-case for AdditionalData::use_global_weights=true and reinit_surface()"));
901
902 do_cell_index_compression =
903 n_unfiltered_cells != numbers::invalid_unsigned_int;
904
905 if (update_flags_mapping & (update_JxW_values | update_normal_vectors))
906 update_flags_mapping |= update_covariant_transformation;
907
908 const unsigned int n_cells = quadrature_vector.size();
909 AssertDimension(n_cells,
910 std::distance(cell_iterator_range.begin(),
911 cell_iterator_range.end()));
912
913 n_q_points_unvectorized.reserve(n_cells);
914
915 // fill unit points index offset vector
916 unit_points_index.reserve(n_cells + 1);
917 unit_points_index.push_back(0);
918 data_index_offsets.reserve(n_cells + 1);
919 data_index_offsets.push_back(0);
920 for (const auto &quadrature : quadrature_vector)
921 {
922 const unsigned int n_points = quadrature.size();
923 n_q_points_unvectorized.push_back(n_points);
924
925 const unsigned int n_q_points =
926 compute_n_q_points<VectorizedArrayType>(n_points);
927 unit_points_index.push_back(unit_points_index.back() + n_q_points);
928
929 const unsigned int n_q_points_data =
930 compute_n_q_points<Number>(n_points);
931 data_index_offsets.push_back(data_index_offsets.back() +
932 n_q_points_data);
933 }
934
935 const unsigned int n_unit_points = unit_points_index.back();
936 const unsigned int n_data_points = data_index_offsets.back();
937
938 // resize data vectors
939 resize_unit_points(n_unit_points);
940 resize_data_fields(n_data_points);
941
942 if (do_cell_index_compression)
943 cell_index_to_compressed_cell_index.resize(n_unfiltered_cells,
945
946 MappingData mapping_data;
947 unsigned int cell_index = 0;
948 for (const auto &cell : cell_iterator_range)
949 {
950 const auto &quadrature = quadrature_vector[cell_index];
951
952 // store unit points
953 const unsigned int n_q_points = compute_n_q_points<VectorizedArrayType>(
954 n_q_points_unvectorized[cell_index]);
955 store_unit_points(unit_points_index[cell_index],
956 n_q_points,
957 n_q_points_unvectorized[cell_index],
958 quadrature_vector[cell_index].get_points());
959
960 // compute mapping data
963 mapping,
964 update_flags_mapping,
965 cell,
966 quadrature,
967 internal_mapping_data,
968 mapping_data);
969
970 // store mapping data
971 const unsigned int n_q_points_data =
972 compute_n_q_points<Number>(n_q_points_unvectorized[cell_index]);
973 store_mapping_data(data_index_offsets[cell_index],
974 n_q_points_data,
975 n_q_points_unvectorized[cell_index],
976 mapping_data,
977 quadrature_vector[cell_index].get_weights());
978
979 if (do_cell_index_compression)
980 cell_index_to_compressed_cell_index[cell->active_cell_index()] =
982
983 ++cell_index;
984 }
985
986 state = State::cell_vector;
987 is_reinitialized();
988 }
989
990
991
992 template <int dim, int spacedim, typename Number>
993 template <typename Iterator>
994 void
996 const IteratorRange<Iterator> & cell_iterator_range,
997 const std::vector<std::vector<Quadrature<dim - 1>>> &quadrature_vector,
998 const unsigned int n_unfiltered_cells)
999 {
1000 clear();
1001
1002 do_cell_index_compression =
1003 n_unfiltered_cells != numbers::invalid_unsigned_int;
1004
1005 const unsigned int n_cells = quadrature_vector.size();
1006 AssertDimension(n_cells,
1007 std::distance(cell_iterator_range.begin(),
1008 cell_iterator_range.end()));
1009
1010 // fill cell index offset vector
1011 cell_index_offset.resize(n_cells);
1012 unsigned int n_faces = 0;
1013 unsigned int cell_index = 0;
1014 for (const auto &cell : cell_iterator_range)
1015 {
1016 cell_index_offset[cell_index] = n_faces;
1017 n_faces += cell->n_faces();
1018 ++cell_index;
1019 }
1020
1021 n_q_points_unvectorized.reserve(n_faces);
1022
1023 // fill unit points index offset vector
1024 unit_points_index.resize(n_faces + 1);
1025 data_index_offsets.resize(n_faces + 1);
1026 cell_index = 0;
1027 unsigned int n_unit_points = 0;
1028 unsigned int n_data_points = 0;
1029 for (const auto &cell : cell_iterator_range)
1030 {
1031 for (const auto &f : cell->face_indices())
1032 {
1033 const unsigned int current_face_index =
1034 cell_index_offset[cell_index] + f;
1035
1036 unit_points_index[current_face_index] = n_unit_points;
1037 data_index_offsets[current_face_index] = n_data_points;
1038
1039 const unsigned int n_points =
1040 quadrature_vector[cell_index][f].size();
1041 n_q_points_unvectorized.push_back(n_points);
1042
1043 const unsigned int n_q_points =
1044 compute_n_q_points<VectorizedArrayType>(n_points);
1045 n_unit_points += n_q_points;
1046
1047 const unsigned int n_q_points_data =
1048 compute_n_q_points<Number>(n_points);
1049 n_data_points += n_q_points_data;
1050 }
1051
1052 ++cell_index;
1053 }
1054 unit_points_index[n_faces] = n_unit_points;
1055 data_index_offsets[n_faces] = n_data_points;
1056
1057 // compress indices
1058 if (do_cell_index_compression)
1059 cell_index_to_compressed_cell_index.resize(n_unfiltered_cells,
1061
1062 // fill unit points and mapping data for every face of all cells
1063 // resize data vectors
1064 resize_unit_points_faces(n_unit_points);
1065 resize_unit_points(n_unit_points);
1066 resize_data_fields(n_data_points);
1067
1068 MappingData mapping_data;
1069 cell_index = 0;
1070 QProjector<dim> q_projector;
1071 for (const auto &cell : cell_iterator_range)
1072 {
1073 const auto &quadratures_on_faces = quadrature_vector[cell_index];
1074
1075 Assert(quadratures_on_faces.size() == cell->n_faces(),
1076 ExcDimensionMismatch(quadratures_on_faces.size(),
1077 cell->n_faces()));
1078
1079 for (const auto &f : cell->face_indices())
1080 {
1081 const auto &quadrature_on_face = quadratures_on_faces[f];
1082
1083 const auto quadrature_on_cell =
1084 q_projector.project_to_face(cell->reference_cell(),
1085 quadrature_on_face,
1086 f);
1087
1088 const unsigned int current_face_index =
1089 cell_index_offset[cell_index] + f;
1090
1091 // store unit points
1092 const unsigned int n_q_points =
1093 compute_n_q_points<VectorizedArrayType>(
1094 n_q_points_unvectorized[current_face_index]);
1095 store_unit_points(unit_points_index[current_face_index],
1096 n_q_points,
1097 n_q_points_unvectorized[current_face_index],
1098 quadrature_on_cell.get_points());
1099
1100 store_unit_points_faces(unit_points_index[current_face_index],
1101 n_q_points,
1102 n_q_points_unvectorized[current_face_index],
1103 quadrature_on_face.get_points());
1104
1107 update_flags_mapping,
1108 cell,
1109 f,
1110 quadrature_on_face,
1111 internal_mapping_data,
1112 mapping_data);
1113
1114 const unsigned int n_q_points_data = compute_n_q_points<Number>(
1115 n_q_points_unvectorized[current_face_index]);
1116 store_mapping_data(data_index_offsets[current_face_index],
1117 n_q_points_data,
1118 n_q_points_unvectorized[current_face_index],
1119 mapping_data,
1120 quadrature_on_face.get_weights());
1121 }
1122 if (do_cell_index_compression)
1123 cell_index_to_compressed_cell_index[cell->active_cell_index()] =
1124 cell_index;
1125
1126 ++cell_index;
1127 }
1128
1129 state = State::faces_on_cells_in_vector;
1130 is_reinitialized();
1131 }
1132
1133
1134
1135 template <int dim, int spacedim, typename Number>
1136 bool
1138 {
1139 return state == State::faces_on_cells_in_vector;
1140 }
1141
1142
1143
1144 template <int dim, int spacedim, typename Number>
1145 unsigned int
1147 const unsigned int cell_index,
1148 const unsigned int face_number) const
1149 {
1151 face_number == numbers::invalid_unsigned_int)
1152 {
1153 Assert(state == State::single_cell,
1154 ExcMessage(
1155 "This mapping info is not reinitialized for a single cell!"));
1156 return n_q_points_unvectorized[0];
1157 }
1158 else
1159 {
1160 return n_q_points_unvectorized[compute_geometry_index_offset(
1161 cell_index, face_number)];
1162 }
1163 }
1164
1165
1166
1167 template <int dim, int spacedim, typename Number>
1168 template <typename NumberType>
1169 unsigned int
1171 const unsigned int n_q_points_unvectorized)
1172 {
1173 const unsigned int n_lanes =
1175 const unsigned int n_filled_lanes_last_batch =
1176 n_q_points_unvectorized % n_lanes;
1177 unsigned int n_q_points = n_q_points_unvectorized / n_lanes;
1178 if (n_filled_lanes_last_batch > 0)
1179 ++n_q_points;
1180 return n_q_points;
1181 }
1182
1183
1184
1185 template <int dim, int spacedim, typename Number>
1186 unsigned int
1188 const unsigned int cell_index,
1189 const unsigned int face_number) const
1190 {
1191 const unsigned int compressed_cell_index =
1192 compute_compressed_cell_index(cell_index);
1193 if (face_number == numbers::invalid_unsigned_int)
1194 {
1195 Assert(state == State::cell_vector,
1196 ExcMessage(
1197 "This mapping info is not reinitialized for a cell vector!"));
1198 return compressed_cell_index;
1199 }
1200 else
1201 {
1203 ExcMessage(
1204 "cell_index has to be set if face_number is specified!"));
1205 Assert(state == State::faces_on_cells_in_vector,
1206 ExcMessage("This mapping info is not reinitialized for faces"
1207 " on cells in a vector!"));
1208 return cell_index_offset[compressed_cell_index] + face_number;
1209 }
1210 }
1211
1212
1213
1214 template <int dim, int spacedim, typename Number>
1215 unsigned int
1217 const unsigned int cell_index) const
1218 {
1219 if (do_cell_index_compression)
1220 {
1221 Assert(cell_index_to_compressed_cell_index[cell_index] !=
1223 ExcMessage("Mapping info object was not initialized for this"
1224 " active cell index!"));
1225 return cell_index_to_compressed_cell_index[cell_index];
1226 }
1227 else
1228 return cell_index;
1229 }
1230
1231
1232 template <int dim, int spacedim, typename Number>
1233 void
1235 const unsigned int unit_points_index_offset,
1236 const unsigned int n_q_points,
1237 const unsigned int n_q_points_unvectorized,
1238 const std::vector<Point<dim>> &points)
1239 {
1240 const unsigned int n_lanes =
1242
1243 for (unsigned int q = 0; q < n_q_points; ++q)
1244 {
1245 const unsigned int offset = unit_points_index_offset + q;
1246 for (unsigned int v = 0;
1247 v < n_lanes && q * n_lanes + v < n_q_points_unvectorized;
1248 ++v)
1249 for (unsigned int d = 0; d < dim; ++d)
1251 unit_points[offset][d], v) = points[q * n_lanes + v][d];
1252 }
1253 }
1254
1255
1256
1257 template <int dim, int spacedim, typename Number>
1258 void
1260 const unsigned int unit_points_index_offset,
1261 const unsigned int n_q_points,
1262 const unsigned int n_q_points_unvectorized,
1263 const std::vector<Point<dim - 1>> &points)
1264 {
1265 const unsigned int n_lanes =
1267
1268 for (unsigned int q = 0; q < n_q_points; ++q)
1269 {
1270 const unsigned int offset = unit_points_index_offset + q;
1271 for (unsigned int v = 0;
1272 v < n_lanes && q * n_lanes + v < n_q_points_unvectorized;
1273 ++v)
1274 for (unsigned int d = 0; d < dim - 1; ++d)
1276 unit_points_faces[offset][d], v) = points[q * n_lanes + v][d];
1277 }
1278 }
1279
1280
1281
1282 template <int dim, int spacedim, typename Number>
1283 void
1285 const unsigned int unit_points_index_offset,
1286 const unsigned int n_q_points,
1287 const unsigned int n_q_points_unvectorized,
1288 const MappingInfo::MappingData &mapping_data,
1289 const std::vector<double> & weights)
1290 {
1291 const unsigned int n_lanes =
1293
1294 for (unsigned int q = 0; q < n_q_points; ++q)
1295 {
1296 const unsigned int offset = unit_points_index_offset + q;
1297 for (unsigned int v = 0;
1298 v < n_lanes && q * n_lanes + v < n_q_points_unvectorized;
1299 ++v)
1300 {
1301 if (update_flags_mapping & UpdateFlags::update_jacobians)
1302 for (unsigned int d = 0; d < dim; ++d)
1303 for (unsigned int s = 0; s < spacedim; ++s)
1305 jacobians[offset][d][s], v) =
1306 mapping_data.jacobians[q * n_lanes + v][d][s];
1307 if (update_flags_mapping & UpdateFlags::update_inverse_jacobians)
1308 for (unsigned int d = 0; d < dim; ++d)
1309 for (unsigned int s = 0; s < spacedim; ++s)
1311 inverse_jacobians[offset][s][d], v) =
1312 mapping_data.inverse_jacobians[q * n_lanes + v][s][d];
1313 if (update_flags_mapping & UpdateFlags::update_JxW_values)
1314 {
1315 if (additional_data.use_global_weights)
1316 {
1318 JxW_values[offset], v) = weights[q * n_lanes + v];
1319 }
1320 else
1321 {
1323 JxW_values[offset], v) =
1324 mapping_data.JxW_values[q * n_lanes + v];
1325 }
1326 }
1327 if (update_flags_mapping & UpdateFlags::update_normal_vectors)
1328 for (unsigned int s = 0; s < spacedim; ++s)
1330 normal_vectors[offset][s], v) =
1331 mapping_data.normal_vectors[q * n_lanes + v][s];
1332 if (update_flags_mapping & UpdateFlags::update_quadrature_points)
1333 for (unsigned int s = 0; s < spacedim; ++s)
1335 real_points[offset][s], v) =
1336 mapping_data.quadrature_points[q * n_lanes + v][s];
1337 }
1338 }
1339 }
1340
1341
1342
1343 template <int dim, int spacedim, typename Number>
1344 void
1346 const unsigned int n_unit_point_batches)
1347 {
1348 unit_points.resize(n_unit_point_batches);
1349 }
1350
1351
1352
1353 template <int dim, int spacedim, typename Number>
1354 void
1356 const unsigned int n_unit_point_batches)
1357 {
1358 unit_points_faces.resize(n_unit_point_batches);
1359 }
1360
1361
1362
1363 template <int dim, int spacedim, typename Number>
1364 void
1366 const unsigned int n_data_point_batches)
1367 {
1368 if (update_flags_mapping & UpdateFlags::update_jacobians)
1369 jacobians.resize(n_data_point_batches);
1370 if (update_flags_mapping & UpdateFlags::update_inverse_jacobians)
1371 inverse_jacobians.resize(n_data_point_batches);
1372 if (update_flags_mapping & UpdateFlags::update_JxW_values)
1373 JxW_values.resize(n_data_point_batches);
1374 if (update_flags_mapping & UpdateFlags::update_normal_vectors)
1375 normal_vectors.resize(n_data_point_batches);
1376 if (update_flags_mapping & UpdateFlags::update_quadrature_points)
1377 real_points.resize(n_data_point_batches);
1378 }
1379
1380
1381
1382 template <int dim, int spacedim, typename Number>
1383 inline const Point<
1384 dim,
1387 const unsigned int offset) const
1388 {
1389 return &unit_points[offset];
1390 }
1391
1392
1393
1394 template <int dim, int spacedim, typename Number>
1395 inline const Point<
1396 dim - 1,
1399 const unsigned int offset) const
1400 {
1401 return &unit_points_faces[offset];
1402 }
1403
1404
1405
1406 template <int dim, int spacedim, typename Number>
1407 inline const Point<dim, Number> *
1409 const unsigned int offset) const
1410 {
1411 return &real_points[offset];
1412 }
1413
1414
1415
1416 template <int dim, int spacedim, typename Number>
1417 unsigned int
1419 const unsigned int cell_index,
1420 const unsigned int face_number) const
1421 {
1423 face_number == numbers::invalid_unsigned_int)
1424 {
1425 Assert(state == State::single_cell,
1426 ExcMessage(
1427 "This mapping info is not reinitialized for a single cell!"));
1428 return 0;
1429 }
1430 else
1431 {
1432 const unsigned int offset =
1433 compute_geometry_index_offset(cell_index, face_number);
1434 return unit_points_index[offset];
1435 }
1436 }
1437
1438
1439
1440 template <int dim, int spacedim, typename Number>
1441 unsigned int
1443 const unsigned int cell_index,
1444 const unsigned int face_number) const
1445 {
1447 face_number == numbers::invalid_unsigned_int)
1448 {
1449 Assert(state == State::single_cell,
1450 ExcMessage(
1451 "This mapping info is not reinitialized for a single cell!"));
1452 return 0;
1453 }
1454 else
1455 {
1456 const unsigned int offset =
1457 compute_geometry_index_offset(cell_index, face_number);
1458 return data_index_offsets[offset];
1459 }
1460 }
1461
1462
1463
1464 template <int dim, int spacedim, typename Number>
1467 const unsigned int offset) const
1468 {
1469 return &jacobians[offset];
1470 }
1471
1472
1473
1474 template <int dim, int spacedim, typename Number>
1477 const unsigned int offset) const
1478 {
1479 return &inverse_jacobians[offset];
1480 }
1481
1482
1483 template <int dim, int spacedim, typename Number>
1484 inline const Tensor<1, spacedim, Number> *
1486 const unsigned int offset) const
1487 {
1488 return &normal_vectors[offset];
1489 }
1490
1491
1492
1493 template <int dim, int spacedim, typename Number>
1494 inline const Number *
1495 MappingInfo<dim, spacedim, Number>::get_JxW(const unsigned int offset) const
1496 {
1497 return &JxW_values[offset];
1498 }
1499
1500
1501
1502 template <int dim, int spacedim, typename Number>
1505 {
1506 return *mapping;
1507 }
1508
1509
1510
1511 template <int dim, int spacedim, typename Number>
1514 {
1515 return update_flags;
1516 }
1517
1518
1519
1520 template <int dim, int spacedim, typename Number>
1523 {
1524 return update_flags_mapping;
1525 }
1526
1527
1528
1529 template <int dim, int spacedim, typename Number>
1530 boost::signals2::connection
1532 const std::function<void()> &set_is_reinitialized)
1533 {
1534 return is_reinitialized.connect(set_is_reinitialized);
1535 }
1536} // namespace NonMatching
1537
1539
1540#endif
IteratorOverIterators end() const
IteratorOverIterators begin()
void initialize_face(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
Definition mapping_q.cc:162
void initialize(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
Definition mapping_q.cc:82
Abstract base class for mapping classes.
Definition mapping.h:317
const UpdateFlags update_flags
AlignedVector< unsigned int > unit_points_index
AlignedVector< Point< dim, VectorizedArrayType > > unit_points
unsigned int compute_unit_point_index_offset(const unsigned int cell_index, const unsigned int face_number) const
const DerivativeForm< 1, spacedim, dim, Number > * get_inverse_jacobian(const unsigned int offset) const
AlignedVector< Number > JxW_values
std::shared_ptr< typename Mapping< dim, spacedim >::InternalDataBase > internal_mapping_data
const Point< dim, Number > * get_real_point(const unsigned int offset) const
boost::signals2::connection connect_is_reinitialized(const std::function< void()> &set_is_reinitialized)
const Mapping< dim, spacedim > & get_mapping() const
AlignedVector< DerivativeForm< 1, dim, spacedim, Number > > jacobians
const Point< dim, VectorizedArrayType > * get_unit_point(const unsigned int offset) const
unsigned int compute_compressed_cell_index(const unsigned int cell_index) const
const Point< dim - 1, VectorizedArrayType > * get_unit_point_faces(const unsigned int offset) const
std::vector< unsigned int > cell_index_offset
void reinit(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Quadrature< dim > &quadrature)
MappingInfo(const Mapping< dim > &mapping, const UpdateFlags update_flags, const AdditionalData additional_data=AdditionalData())
const SmartPointer< const Mapping< dim, spacedim > > mapping
AlignedVector< Tensor< 1, spacedim, Number > > normal_vectors
void reinit(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const ArrayView< const Point< dim > > &unit_points)
void reinit_cells(const ContainerType &cell_iterator_range, const std::vector< Quadrature< dim > > &quadrature_vector, const unsigned int n_unfiltered_cells=numbers::invalid_unsigned_int)
const DerivativeForm< 1, dim, spacedim, Number > * get_jacobian(const unsigned int offset) const
boost::signals2::signal< void()> is_reinitialized
unsigned int get_n_q_points_unvectorized(const unsigned int cell_index, const unsigned int face_number) const
void resize_unit_points_faces(const unsigned int n_unit_point_batches)
unsigned int compute_data_index_offset(const unsigned int cell_index, const unsigned int face_number) const
const Tensor< 1, spacedim, Number > * get_normal_vector(const unsigned int offset) const
void reinit(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const std::vector< Point< dim > > &unit_points)
void store_mapping_data(const unsigned int unit_points_index_offset, const unsigned int n_q_points, const unsigned int n_q_points_unvectorized, const MappingData &mapping_data, const std::vector< double > &weights)
typename ::internal::VectorizedArrayTrait< Number >::vectorized_value_type VectorizedArrayType
void reinit_cells(const ContainerType &cell_iterator_range, const std::vector< std::vector< Point< dim > > > &unit_points_vector, const unsigned int n_unfiltered_cells=numbers::invalid_unsigned_int)
void resize_unit_points(const unsigned int n_unit_point_batches)
UpdateFlags get_update_flags() const
void store_unit_points_faces(const unsigned int unit_points_index_offset, const unsigned int n_q_points, const unsigned int n_q_points_unvectorized, const std::vector< Point< dim - 1 > > &points)
void resize_data_fields(const unsigned int n_data_point_batches)
std::vector< unsigned int > cell_index_to_compressed_cell_index
AlignedVector< Point< spacedim, Number > > real_points
unsigned int compute_geometry_index_offset(const unsigned int cell_index, const unsigned int face_number) const
AlignedVector< DerivativeForm< 1, spacedim, dim, Number > > inverse_jacobians
std::vector< unsigned int > n_q_points_unvectorized
const Number * get_JxW(const unsigned int offset) const
void store_unit_points(const unsigned int unit_points_index_offset, const unsigned int n_q_points, const unsigned int n_q_points_unvectorized, const std::vector< Point< dim > > &points)
AlignedVector< unsigned int > data_index_offsets
const AdditionalData additional_data
AlignedVector< Point< dim - 1, VectorizedArrayType > > unit_points_faces
unsigned int compute_n_q_points(const unsigned int n_q_points_unvectorized)
void reinit_surface(const IteratorRange< Iterator > &cell_iterator_range, const std::vector< ImmersedSurfaceQuadrature< dim > > &quadrature_vector, const unsigned int n_unfiltered_cells=numbers::invalid_unsigned_int)
void reinit_faces(const IteratorRange< Iterator > &cell_iterator_range, const std::vector< std::vector< Quadrature< dim - 1 > > > &quadrature_vector, const unsigned int n_unfiltered_cells=numbers::invalid_unsigned_int)
UpdateFlags update_flags_mapping
UpdateFlags get_update_flags_mapping() const
static void compute_mapping_data_for_immersed_surface_quadrature(const SmartPointer< const Mapping< dim, spacedim > > mapping, const UpdateFlags &update_flags_mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const ImmersedSurfaceQuadrature< dim > &quadrature, std::shared_ptr< typename Mapping< dim, spacedim >::InternalDataBase > internal_mapping_data, MappingData &mapping_data)
static void compute_mapping_data_for_quadrature(const SmartPointer< const Mapping< dim, spacedim > > mapping, const UpdateFlags &update_flags_mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, CellSimilarity::Similarity &cell_similarity, const Quadrature< dim > &quadrature, std::shared_ptr< typename Mapping< dim, spacedim >::InternalDataBase > internal_mapping_data, MappingData &mapping_data)
static void compute_mapping_data_for_face_quadrature(const SmartPointer< const Mapping< dim, spacedim > > mapping, const UpdateFlags &update_flags_mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< dim - 1 > &quadrature, std::shared_ptr< typename Mapping< dim, spacedim >::InternalDataBase > internal_mapping_data, MappingData &mapping_data)
static UpdateFlags required_update_flags(const SmartPointer< const Mapping< dim, spacedim > > mapping, const UpdateFlags &update_flags)
Definition point.h:112
static void project_to_face(const ReferenceCell &reference_cell, const SubQuadrature &quadrature, const unsigned int face_no, std::vector< Point< dim > > &q_points)
const std::vector< double > & get_weights() const
const std::vector< Point< dim > > & get_points() const
unsigned int size() const
std::vector< DerivativeForm< 1, spacedim, dim > > inverse_jacobians
void initialize(const unsigned int n_quadrature_points, const UpdateFlags flags)
std::vector< DerivativeForm< 1, dim, spacedim > > jacobians
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
unsigned int cell_index
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcMessage(std::string arg1)
UpdateFlags
@ 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_gradients
Shape function gradients.
@ update_quadrature_points
Transformed quadrature points.
@ update_default
No update.
static const unsigned int invalid_unsigned_int
Definition types.h:213
AdditionalData(const bool use_global_weights=false)
static constexpr std::size_t width()
static value_type & get(value_type &value, unsigned int c)