Reference documentation for deal.II version GIT d6cf33b98c 2023-09-22 19:50:02+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\}}\)
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 
27 
28 #include <deal.II/fe/fe_dgq.h>
30 #include <deal.II/fe/mapping.h>
32 #include <deal.II/fe/mapping_q.h>
34 
36 
37 #include <memory>
38 
39 
41 
42 namespace NonMatching
43 {
44  namespace internal
45  {
46  template <int dim, int spacedim = dim>
48  {
49  using MappingData =
51  spacedim>;
52 
53  public:
54  static UpdateFlags
56  const SmartPointer<const Mapping<dim, spacedim>> mapping,
57  const UpdateFlags &update_flags)
58  {
59  return mapping->requires_update_flags(update_flags);
60  }
61 
62  static void
64  const SmartPointer<const Mapping<dim, spacedim>> mapping,
65  const UpdateFlags &update_flags_mapping,
67  CellSimilarity::Similarity &cell_similarity,
68  const Quadrature<dim> &quadrature,
69  std::shared_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
70  internal_mapping_data,
71  MappingData &mapping_data)
72  {
73  mapping_data.initialize(quadrature.size(), update_flags_mapping);
74 
75  // reuse internal_mapping_data for MappingQ to avoid memory allocations
76  if (const MappingQ<dim, spacedim> *mapping_q =
77  dynamic_cast<const MappingQ<dim, spacedim> *>(&(*mapping)))
78  {
79  (void)mapping_q;
80  auto &data =
81  dynamic_cast<typename MappingQ<dim, spacedim>::InternalData &>(
82  *internal_mapping_data);
83  data.initialize(update_flags_mapping,
84  quadrature,
85  quadrature.size());
86  }
87  else
88  {
89  internal_mapping_data =
90  mapping->get_data(update_flags_mapping, quadrature);
91  }
92 
93  cell_similarity = mapping->fill_fe_values(cell,
94  cell_similarity,
95  quadrature,
96  *internal_mapping_data,
97  mapping_data);
98  }
99 
100 
101 
102  static void
104  const SmartPointer<const Mapping<dim, spacedim>> mapping,
105  const UpdateFlags &update_flags_mapping,
106  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
107  const ImmersedSurfaceQuadrature<dim> &quadrature,
108  std::shared_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
109  internal_mapping_data,
110  MappingData &mapping_data)
111  {
112  mapping_data.initialize(quadrature.size(), update_flags_mapping);
113 
114  // reuse internal_mapping_data for MappingQ to avoid memory allocations
115  if (dynamic_cast<const MappingQ<dim, spacedim> *>(&(*mapping)))
116  {
117  auto &data =
118  dynamic_cast<typename MappingQ<dim, spacedim>::InternalData &>(
119  *internal_mapping_data);
120  data.initialize(update_flags_mapping,
121  quadrature,
122  quadrature.size());
123  }
124  else
125  {
126  internal_mapping_data =
127  mapping->get_data(update_flags_mapping, quadrature);
128  }
129 
130  mapping->fill_fe_immersed_surface_values(cell,
131  quadrature,
132  *internal_mapping_data,
133  mapping_data);
134  }
135 
136 
137 
138  static void
140  const SmartPointer<const Mapping<dim, spacedim>> mapping,
141  const UpdateFlags &update_flags_mapping,
142  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
143  const unsigned int face_no,
144  const Quadrature<dim - 1> &quadrature,
145  std::shared_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
146  internal_mapping_data,
147  MappingData &mapping_data)
148  {
149  mapping_data.initialize(quadrature.size(), update_flags_mapping);
150 
151  // reuse internal_mapping_data for MappingQ to avoid memory allocations
152  if (const MappingQ<dim, spacedim> *mapping_q =
153  dynamic_cast<const MappingQ<dim, spacedim> *>(&(*mapping)))
154  {
155  auto &data =
156  dynamic_cast<typename MappingQ<dim, spacedim>::InternalData &>(
157  *internal_mapping_data);
158  data.initialize_face(update_flags_mapping,
160  ReferenceCells::get_hypercube<dim>(),
161  quadrature,
162  face_no,
163  cell->face_orientation(face_no),
164  cell->face_flip(face_no),
165  cell->face_rotation(face_no)),
166  quadrature.size());
167 
168  mapping_q->fill_mapping_data_for_face_quadrature(
169  cell, face_no, quadrature, *internal_mapping_data, mapping_data);
170  }
171  else
172  {
173  auto internal_mapping_data =
174  mapping->get_face_data(update_flags_mapping,
175  hp::QCollection<dim - 1>(quadrature));
176 
177  mapping->fill_fe_face_values(cell,
178  face_no,
179  hp::QCollection<dim - 1>(quadrature),
180  *internal_mapping_data,
181  mapping_data);
182  }
183  }
184  };
185 
186  template <int dim, int spacedim = dim>
189  const double diameter,
191  &inverse_jacobians)
192  {
193  const auto jac_0 = inverse_jacobians[0];
194  const double zero_tolerance_double =
196  bool jacobian_constant = true;
197  for (unsigned int q = 1; q < inverse_jacobians.size(); ++q)
198  {
199  const DerivativeForm<1, dim, spacedim> &jac = inverse_jacobians[q];
200  for (unsigned int d = 0; d < dim; ++d)
201  for (unsigned int e = 0; e < spacedim; ++e)
202  if (std::fabs(jac_0[d][e] - jac[d][e]) > zero_tolerance_double)
203  jacobian_constant = false;
204  if (!jacobian_constant)
205  break;
206  }
207 
208  // check whether the Jacobian is diagonal to machine
209  // accuracy
210  bool cell_cartesian = jacobian_constant;
211  for (unsigned int d = 0; d < dim; ++d)
212  for (unsigned int e = 0; e < dim; ++e)
213  if (d != e)
214  if (std::fabs(jac_0[d][e]) > zero_tolerance_double)
215  {
216  cell_cartesian = false;
217  break;
218  }
219 
220  // return cell type
221  if (cell_cartesian)
223  else if (jacobian_constant)
225  else
227  }
228  } // namespace internal
229 
246  template <int dim, int spacedim = dim, typename Number = double>
247  class MappingInfo : public Subscriptor
248  {
249  public:
253  using VectorizedArrayType = typename ::internal::VectorizedArrayTrait<
254  Number>::vectorized_value_type;
255 
261  {
265  AdditionalData(const bool use_global_weights = false,
266  const bool store_cells = false)
269  {}
270 
277 
286  };
287 
303 
308  void
310  const std::vector<Point<dim>> &unit_points);
311 
316  void
318  const ArrayView<const Point<dim>> &unit_points);
319 
325  void
327  const Quadrature<dim> &quadrature);
328 
339  template <typename ContainerType>
340  void
341  reinit_cells(
342  const ContainerType &cell_iterator_range,
343  const std::vector<std::vector<Point<dim>>> &unit_points_vector,
344  const unsigned int n_unfiltered_cells = numbers::invalid_unsigned_int);
345 
352  template <typename ContainerType>
353  void
354  reinit_cells(
355  const ContainerType &cell_iterator_range,
356  const std::vector<Quadrature<dim>> &quadrature_vector,
357  const unsigned int n_unfiltered_cells = numbers::invalid_unsigned_int);
358 
363  template <typename Iterator>
364  void
366  const IteratorRange<Iterator> &cell_iterator_range,
367  const std::vector<ImmersedSurfaceQuadrature<dim>> &quadrature_vector,
368  const unsigned int n_unfiltered_cells = numbers::invalid_unsigned_int);
369 
374  template <typename Iterator>
375  void
376  reinit_faces(
377  const IteratorRange<Iterator> &cell_iterator_range,
378  const std::vector<std::vector<Quadrature<dim - 1>>> &quadrature_vector,
379  const unsigned int n_unfiltered_cells = numbers::invalid_unsigned_int);
380 
385  bool
386  is_face_state() const;
387 
393  get_unit_point(const unsigned int offset) const;
394 
399  const Point<dim - 1, VectorizedArrayType> *
400  get_unit_point_faces(const unsigned int offset) const;
401 
407  get_jacobian(const unsigned int offset) const;
408 
414  get_inverse_jacobian(const unsigned int offset) const;
415 
421  get_normal_vector(const unsigned int offset) const;
422 
427  const Number *
428  get_JxW(const unsigned int offset) const;
429 
434  const Point<dim, Number> *
435  get_real_point(const unsigned int offset) const;
436 
440  const Mapping<dim, spacedim> &
441  get_mapping() const;
442 
447  get_update_flags() const;
448 
453  get_update_flags_mapping() const;
454 
458  boost::signals2::connection
459  connect_is_reinitialized(const std::function<void()> &set_is_reinitialized);
460 
464  template <bool is_face>
465  unsigned int
466  compute_geometry_index_offset(const unsigned int cell_index,
467  const unsigned int face_number) const;
468 
472  unsigned int
473  compute_unit_point_index_offset(const unsigned int geometry_index) const;
474 
478  unsigned int
479  compute_data_index_offset(const unsigned int geometry_index) const;
480 
484  unsigned int
486  const unsigned int geometry_index) const;
487 
491  unsigned int
492  get_n_q_points_unvectorized(const unsigned int geometry_index) const;
493 
498  get_cell_type(const unsigned int geometry_index) const;
499 
506  get_cell_iterator(const unsigned int cell_index) const;
507 
511  std::size_t
512  memory_consumption() const;
513 
514  private:
515  using MappingData =
517  spacedim>;
518 
522  template <typename NumberType>
523  unsigned int
524  compute_n_q_points(const unsigned int n_q_points_unvectorized);
525 
529  void
530  clear();
531 
535  void
536  resize_unit_points(const unsigned int n_unit_point_batches);
537 
541  void
542  resize_unit_points_faces(const unsigned int n_unit_point_batches);
543 
547  void
548  resize_data_fields(const unsigned int n_data_point_batches);
549 
553  void
554  store_unit_points(const unsigned int unit_points_index_offset,
555  const unsigned int n_q_points,
556  const unsigned int n_q_points_unvectorized,
557  const std::vector<Point<dim>> &points);
558 
562  void
563  store_unit_points_faces(const unsigned int unit_points_index_offset,
564  const unsigned int n_q_points,
565  const unsigned int n_q_points_unvectorized,
566  const std::vector<Point<dim - 1>> &points);
567 
571  void
572  store_mapping_data(const unsigned int unit_points_index_offset,
573  const unsigned int n_q_points,
574  const unsigned int n_q_points_unvectorized,
575  const MappingData &mapping_data,
576  const std::vector<double> &weights,
577  const unsigned int compressed_unit_point_index_offset,
578  const bool affine_cell);
579 
583  unsigned int
584  compute_compressed_cell_index(const unsigned int cell_index) const;
585 
589  template <typename ContainerType, typename QuadratureType>
590  void
592  const ContainerType &cell_iterator_range,
593  const std::vector<QuadratureType> &quadrature_vector,
594  const unsigned int n_unfiltered_cells,
595  const std::function<
596  void(const typename Triangulation<dim, spacedim>::cell_iterator &cell,
597  const QuadratureType &quadrature,
598  MappingData &mapping_data)> &compute_mapping_data);
599 
603  enum class State
604  {
605  invalid,
606  single_cell,
607  cell_vector,
609  };
610 
616 
623 
630 
634  std::vector<unsigned int> unit_points_index;
635 
639  std::shared_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
641 
646 
651 
656 
661 
668  std::vector<::internal::MatrixFreeFunctions::GeometryType> cell_type;
669 
673  std::vector<unsigned int> data_index_offsets;
674 
678  std::vector<unsigned int> compressed_data_index_offsets;
679 
687 
694 
702 
711 
718 
722  std::vector<unsigned int> n_q_points_unvectorized;
723 
728  std::vector<unsigned int> cell_index_offset;
729 
734  std::vector<unsigned int> cell_index_to_compressed_cell_index;
735 
740 
745  boost::signals2::signal<void()> is_reinitialized;
746 
753 
758  std::vector<std::pair<int, int>> cell_level_and_indices;
759  };
760 
761  // ----------------------- template functions ----------------------
762 
763 
764  template <int dim, int spacedim, typename Number>
766  const Mapping<dim> &mapping,
767  const UpdateFlags update_flags,
768  const AdditionalData additional_data)
769  : mapping(&mapping)
770  , update_flags(update_flags)
771  , update_flags_mapping(update_default)
772  , additional_data(additional_data)
773  {
774  // translate update flags
784 
785  // always save quadrature points for now
787 
790  this->mapping, update_flags_mapping);
791 
792  // construct internal_mapping_data for MappingQ to be able to reuse it in
793  // reinit() calls to avoid memory allocations
794  if (const MappingQ<dim, spacedim> *mapping_q =
795  dynamic_cast<const MappingQ<dim, spacedim> *>(&mapping))
796  {
798  std::make_unique<typename MappingQ<dim, spacedim>::InternalData>(
799  mapping_q->get_degree());
800  }
801  }
802 
803 
804 
805  template <int dim, int spacedim, typename Number>
806  void
808  {
809  n_q_points_unvectorized.clear();
810  unit_points_index.clear();
811  data_index_offsets.clear();
812  compressed_data_index_offsets.clear();
813  cell_type.clear();
814  }
815 
816 
817 
818  template <int dim, int spacedim, typename Number>
819  void
821  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
822  const std::vector<Point<dim>> &unit_points_in)
823  {
824  reinit(cell, Quadrature<dim>(unit_points_in));
825  }
826 
827 
828 
829  template <int dim, int spacedim, typename Number>
830  void
832  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
833  const ArrayView<const Point<dim>> &unit_points_in)
834  {
835  reinit(cell,
836  std::vector<Point<dim>>(unit_points_in.begin(),
837  unit_points_in.end()));
838  }
839 
840 
841 
842  template <int dim, int spacedim, typename Number>
843  void
845  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
846  const Quadrature<dim> &quadrature)
847  {
848  n_q_points_unvectorized.resize(1);
849  n_q_points_unvectorized[0] = quadrature.size();
850 
851  const unsigned int n_q_points =
852  compute_n_q_points<VectorizedArrayType>(n_q_points_unvectorized[0]);
853 
854  const unsigned int n_q_points_data =
855  compute_n_q_points<Number>(n_q_points_unvectorized[0]);
856 
857  // resize data vectors
858  resize_unit_points(n_q_points);
859  resize_data_fields(n_q_points_data);
860 
861  // store unit points
862  store_unit_points(0,
863  n_q_points,
864  n_q_points_unvectorized[0],
865  quadrature.get_points());
866 
867  // compute mapping data
868  MappingData mapping_data;
872  update_flags_mapping,
873  cell,
874  cell_similarity,
875  quadrature,
876  internal_mapping_data,
877  mapping_data);
878 
879  // check for cartesian/affine cell
880  if (!quadrature.empty() &&
881  update_flags_mapping & UpdateFlags::update_inverse_jacobians)
882  {
883  cell_type.push_back(
884  internal::compute_geometry_type(cell->diameter(),
885  mapping_data.inverse_jacobians));
886  }
887  else
888  cell_type.push_back(
890 
891  // store mapping data
892  store_mapping_data(
893  0,
894  n_q_points_data,
895  n_q_points_unvectorized[0],
896  mapping_data,
897  quadrature.get_weights(),
898  0,
899  cell_type.back() <=
901 
902  unit_points_index.push_back(0);
903  data_index_offsets.push_back(0);
904  compressed_data_index_offsets.push_back(0);
905 
906  state = State::single_cell;
907  is_reinitialized();
908  }
909 
910 
911 
912  template <int dim, int spacedim, typename Number>
913  template <typename ContainerType>
914  void
916  const ContainerType &cell_iterator_range,
917  const std::vector<std::vector<Point<dim>>> &unit_points_vector,
918  const unsigned int n_unfiltered_cells)
919  {
920  const unsigned int n_cells = unit_points_vector.size();
922  std::distance(cell_iterator_range.begin(),
923  cell_iterator_range.end()));
924 
925  std::vector<Quadrature<dim>> quadrature_vector(n_cells);
926  for (unsigned int cell_index = 0; cell_index < n_cells; ++cell_index)
927  quadrature_vector[cell_index] =
928  Quadrature<dim>(unit_points_vector[cell_index]);
929 
930  reinit_cells(cell_iterator_range, quadrature_vector, n_unfiltered_cells);
931  }
932 
933 
934 
935  template <int dim, int spacedim, typename Number>
936  template <typename ContainerType, typename QuadratureType>
937  void
939  const ContainerType &cell_iterator_range,
940  const std::vector<QuadratureType> &quadrature_vector,
941  const unsigned int n_unfiltered_cells,
942  const std::function<
943  void(const typename Triangulation<dim, spacedim>::cell_iterator &cell,
944  const QuadratureType &quadrature,
945  MappingData &mapping_data)> &compute_mapping_data)
946  {
947  clear();
948 
949  do_cell_index_compression =
950  n_unfiltered_cells != numbers::invalid_unsigned_int;
951 
952  const unsigned int n_cells = quadrature_vector.size();
954  std::distance(cell_iterator_range.begin(),
955  cell_iterator_range.end()));
956 
957  n_q_points_unvectorized.reserve(n_cells);
958 
959  cell_type.reserve(n_cells);
960 
961  if (additional_data.store_cells)
962  cell_level_and_indices.resize(n_cells);
963 
964  // fill unit points index offset vector
965  unit_points_index.reserve(n_cells + 1);
966  unit_points_index.push_back(0);
967  data_index_offsets.reserve(n_cells + 1);
968  data_index_offsets.push_back(0);
969  for (const auto &quadrature : quadrature_vector)
970  {
971  const unsigned int n_points = quadrature.size();
972  n_q_points_unvectorized.push_back(n_points);
973 
974  const unsigned int n_q_points =
975  compute_n_q_points<VectorizedArrayType>(n_points);
976  unit_points_index.push_back(unit_points_index.back() + n_q_points);
977 
978  const unsigned int n_q_points_data =
979  compute_n_q_points<Number>(n_points);
980  data_index_offsets.push_back(data_index_offsets.back() +
981  n_q_points_data);
982  }
983 
984  const unsigned int n_unit_points = unit_points_index.back();
985  const unsigned int n_data_points = data_index_offsets.back();
986 
987  // resize data vectors
988  resize_unit_points(n_unit_points);
989  resize_data_fields(n_data_points);
990 
991  if (do_cell_index_compression)
992  cell_index_to_compressed_cell_index.resize(n_unfiltered_cells,
994 
995  MappingData mapping_data;
996  MappingData mapping_data_last_cell;
997  unsigned int size_compressed_data = 0;
998  unsigned int cell_index = 0;
999  for (const auto &cell : cell_iterator_range)
1000  {
1001  if (additional_data.store_cells)
1002  {
1003  this->triangulation = &cell->get_triangulation();
1004  cell_level_and_indices[cell_index] = {cell->level(), cell->index()};
1005  }
1006 
1007  const auto &quadrature = quadrature_vector[cell_index];
1008  const bool empty = quadrature.empty();
1009 
1010  // store unit points
1011  const unsigned int n_q_points = compute_n_q_points<VectorizedArrayType>(
1012  n_q_points_unvectorized[cell_index]);
1013  store_unit_points(unit_points_index[cell_index],
1014  n_q_points,
1015  n_q_points_unvectorized[cell_index],
1016  quadrature.get_points());
1017 
1018  // compute mapping data
1019  compute_mapping_data(cell, quadrature, mapping_data);
1020 
1021  // store mapping data
1022  const unsigned int n_q_points_data =
1023  compute_n_q_points<Number>(n_q_points_unvectorized[cell_index]);
1024 
1025  // check for cartesian/affine cell
1026  if (!empty &&
1027  update_flags_mapping & UpdateFlags::update_inverse_jacobians)
1028  {
1029  cell_type.push_back(
1030  internal::compute_geometry_type(cell->diameter(),
1031  mapping_data.inverse_jacobians));
1032  }
1033  else
1034  cell_type.push_back(
1036 
1037  if (cell_index > 0)
1038  {
1039  // check if current and last cell are affine
1040  const bool affine_cells =
1041  cell_type[cell_index] <=
1043  cell_type[cell_index - 1] <=
1045 
1046  // create a comparator to compare inverse Jacobian of current
1047  // and last cell
1049  1e4 / cell->diameter() * std::numeric_limits<double>::epsilon() *
1050  1024.);
1051 
1052  // we can only compare if current and last cell have at least
1053  // one quadrature point and both cells are at least affine
1054  const auto comparison_result =
1055  (!affine_cells || mapping_data.inverse_jacobians.empty() ||
1056  mapping_data_last_cell.inverse_jacobians.empty()) ?
1058  comparator.compare(mapping_data.inverse_jacobians[0],
1059  mapping_data_last_cell.inverse_jacobians[0]);
1060 
1061  // we can compress the Jacobians and inverse Jacobians if
1062  // inverse Jacobians are equal and cells are affine
1063  if (affine_cells &&
1064  comparison_result ==
1066  {
1067  compressed_data_index_offsets.push_back(
1068  compressed_data_index_offsets.back());
1069  }
1070  else
1071  {
1072  const unsigned int n_compressed_data_last_cell =
1073  cell_type[cell_index - 1] <=
1075  1 :
1076  compute_n_q_points<Number>(
1077  n_q_points_unvectorized[cell_index - 1]);
1078 
1079  compressed_data_index_offsets.push_back(
1080  compressed_data_index_offsets.back() +
1081  n_compressed_data_last_cell);
1082  }
1083  }
1084  else
1085  compressed_data_index_offsets.push_back(0);
1086 
1087  // cache mapping_data from last cell
1088  mapping_data_last_cell = mapping_data;
1089 
1090  store_mapping_data(data_index_offsets[cell_index],
1091  n_q_points_data,
1092  n_q_points_unvectorized[cell_index],
1093  mapping_data,
1094  quadrature.get_weights(),
1095  compressed_data_index_offsets[cell_index],
1096  cell_type[cell_index] <=
1098 
1099  // update size of compressed data depending on cell type and handle
1100  // empty quadratures
1101  if (cell_type[cell_index] <=
1103  size_compressed_data = compressed_data_index_offsets.back() + 1;
1104  else
1105  size_compressed_data =
1106  std::max(size_compressed_data,
1107  compressed_data_index_offsets.back() + n_q_points_data);
1108 
1109  if (do_cell_index_compression)
1110  cell_index_to_compressed_cell_index[cell->active_cell_index()] =
1111  cell_index;
1112 
1113  ++cell_index;
1114  }
1115 
1116  if (update_flags_mapping & UpdateFlags::update_jacobians)
1117  {
1118  jacobians.resize(size_compressed_data);
1119  jacobians.shrink_to_fit();
1120  }
1121  if (update_flags_mapping & UpdateFlags::update_inverse_jacobians)
1122  {
1123  inverse_jacobians.resize(size_compressed_data);
1124  inverse_jacobians.shrink_to_fit();
1125  }
1126 
1127  state = State::cell_vector;
1128  is_reinitialized();
1129  }
1130 
1131 
1132 
1133  template <int dim, int spacedim, typename Number>
1134  template <typename ContainerType>
1135  void
1137  const ContainerType &cell_iterator_range,
1138  const std::vector<Quadrature<dim>> &quadrature_vector,
1139  const unsigned int n_unfiltered_cells)
1140  {
1141  auto compute_mapping_data_for_cells =
1142  [&](const typename Triangulation<dim, spacedim>::cell_iterator &cell,
1143  const Quadrature<dim> &quadrature,
1144  MappingData &mapping_data) {
1145  CellSimilarity::Similarity cell_similarity =
1149  update_flags_mapping,
1150  cell,
1151  cell_similarity,
1152  quadrature,
1153  internal_mapping_data,
1154  mapping_data);
1155  };
1156 
1157  do_reinit_cells<ContainerType, Quadrature<dim>>(
1158  cell_iterator_range,
1159  quadrature_vector,
1160  n_unfiltered_cells,
1161  compute_mapping_data_for_cells);
1162  }
1163 
1164 
1165 
1166  template <int dim, int spacedim, typename Number>
1167  template <typename Iterator>
1168  void
1170  const IteratorRange<Iterator> &cell_iterator_range,
1171  const std::vector<ImmersedSurfaceQuadrature<dim>> &quadrature_vector,
1172  const unsigned int n_unfiltered_cells)
1173  {
1174  Assert(
1175  additional_data.use_global_weights == false,
1176  ExcMessage(
1177  "There is no known use-case for AdditionalData::use_global_weights=true and reinit_surface()"));
1178 
1179  Assert(additional_data.store_cells == false, ExcNotImplemented());
1180 
1181  if (update_flags_mapping & (update_JxW_values | update_normal_vectors))
1182  update_flags_mapping |= update_covariant_transformation;
1183 
1184  auto compute_mapping_data_for_surface =
1185  [&](const typename Triangulation<dim, spacedim>::cell_iterator &cell,
1186  const ImmersedSurfaceQuadrature<dim> &quadrature,
1187  MappingData &mapping_data) {
1190  mapping,
1191  update_flags_mapping,
1192  cell,
1193  quadrature,
1194  internal_mapping_data,
1195  mapping_data);
1196  };
1197 
1198  do_reinit_cells<IteratorRange<Iterator>, ImmersedSurfaceQuadrature<dim>>(
1199  cell_iterator_range,
1200  quadrature_vector,
1201  n_unfiltered_cells,
1202  compute_mapping_data_for_surface);
1203  }
1204 
1205 
1206 
1207  template <int dim, int spacedim, typename Number>
1208  template <typename Iterator>
1209  void
1211  const IteratorRange<Iterator> &cell_iterator_range,
1212  const std::vector<std::vector<Quadrature<dim - 1>>> &quadrature_vector,
1213  const unsigned int n_unfiltered_cells)
1214  {
1215  clear();
1216 
1217  Assert(additional_data.store_cells == false, ExcNotImplemented());
1218 
1219  do_cell_index_compression =
1220  n_unfiltered_cells != numbers::invalid_unsigned_int;
1221 
1222  const unsigned int n_cells = quadrature_vector.size();
1224  std::distance(cell_iterator_range.begin(),
1225  cell_iterator_range.end()));
1226 
1227  // fill cell index offset vector
1228  cell_index_offset.resize(n_cells);
1229  unsigned int n_faces = 0;
1230  unsigned int cell_index = 0;
1231  for (const auto &cell : cell_iterator_range)
1232  {
1233  cell_index_offset[cell_index] = n_faces;
1234  n_faces += cell->n_faces();
1235  ++cell_index;
1236  }
1237 
1238  n_q_points_unvectorized.reserve(n_faces);
1239 
1240  cell_type.reserve(n_faces);
1241 
1242  // fill unit points index offset vector
1243  unit_points_index.resize(n_faces + 1);
1244  data_index_offsets.resize(n_faces + 1);
1245  cell_index = 0;
1246  unsigned int n_unit_points = 0;
1247  unsigned int n_data_points = 0;
1248  for (const auto &cell : cell_iterator_range)
1249  {
1250  for (const auto &f : cell->face_indices())
1251  {
1252  const unsigned int current_face_index =
1253  cell_index_offset[cell_index] + f;
1254 
1255  unit_points_index[current_face_index] = n_unit_points;
1256  data_index_offsets[current_face_index] = n_data_points;
1257 
1258  const unsigned int n_points =
1259  quadrature_vector[cell_index][f].size();
1260  n_q_points_unvectorized.push_back(n_points);
1261 
1262  const unsigned int n_q_points =
1263  compute_n_q_points<VectorizedArrayType>(n_points);
1264  n_unit_points += n_q_points;
1265 
1266  const unsigned int n_q_points_data =
1267  compute_n_q_points<Number>(n_points);
1268  n_data_points += n_q_points_data;
1269  }
1270 
1271  ++cell_index;
1272  }
1273  unit_points_index[n_faces] = n_unit_points;
1274  data_index_offsets[n_faces] = n_data_points;
1275 
1276  // compress indices
1277  if (do_cell_index_compression)
1278  cell_index_to_compressed_cell_index.resize(n_unfiltered_cells,
1280 
1281  // fill unit points and mapping data for every face of all cells
1282  // resize data vectors
1283  resize_unit_points_faces(n_unit_points);
1284  resize_unit_points(n_unit_points);
1285  resize_data_fields(n_data_points);
1286 
1287  MappingData mapping_data;
1288  cell_index = 0;
1289  QProjector<dim> q_projector;
1290  for (const auto &cell : cell_iterator_range)
1291  {
1292  const auto &quadratures_on_faces = quadrature_vector[cell_index];
1293 
1294  Assert(quadratures_on_faces.size() == cell->n_faces(),
1295  ExcDimensionMismatch(quadratures_on_faces.size(),
1296  cell->n_faces()));
1297 
1298  for (const auto &f : cell->face_indices())
1299  {
1300  const auto &quadrature_on_face = quadratures_on_faces[f];
1301 
1302  const auto quadrature_on_cell =
1303  q_projector.project_to_face(cell->reference_cell(),
1304  quadrature_on_face,
1305  f);
1306 
1307  const unsigned int current_face_index =
1308  cell_index_offset[cell_index] + f;
1309 
1310  // store unit points
1311  const unsigned int n_q_points =
1312  compute_n_q_points<VectorizedArrayType>(
1313  n_q_points_unvectorized[current_face_index]);
1314  store_unit_points(unit_points_index[current_face_index],
1315  n_q_points,
1316  n_q_points_unvectorized[current_face_index],
1317  quadrature_on_cell.get_points());
1318 
1319  store_unit_points_faces(unit_points_index[current_face_index],
1320  n_q_points,
1321  n_q_points_unvectorized[current_face_index],
1322  quadrature_on_face.get_points());
1323 
1326  update_flags_mapping,
1327  cell,
1328  f,
1329  quadrature_on_face,
1330  internal_mapping_data,
1331  mapping_data);
1332 
1333  cell_type.push_back(
1335 
1336  compressed_data_index_offsets.push_back(
1337  data_index_offsets[current_face_index]);
1338 
1339  const unsigned int n_q_points_data = compute_n_q_points<Number>(
1340  n_q_points_unvectorized[current_face_index]);
1341  store_mapping_data(data_index_offsets[current_face_index],
1342  n_q_points_data,
1343  n_q_points_unvectorized[current_face_index],
1344  mapping_data,
1345  quadrature_on_face.get_weights(),
1346  data_index_offsets[current_face_index],
1347  false);
1348  }
1349  if (do_cell_index_compression)
1350  cell_index_to_compressed_cell_index[cell->active_cell_index()] =
1351  cell_index;
1352 
1353  ++cell_index;
1354  }
1355 
1356  state = State::faces_on_cells_in_vector;
1357  is_reinitialized();
1358  }
1359 
1360 
1361 
1362  template <int dim, int spacedim, typename Number>
1363  bool
1365  {
1366  return state == State::faces_on_cells_in_vector;
1367  }
1368 
1369 
1370 
1371  template <int dim, int spacedim, typename Number>
1372  unsigned int
1374  const unsigned int geometry_index) const
1375  {
1376  return n_q_points_unvectorized[geometry_index];
1377  }
1378 
1379 
1380  template <int dim, int spacedim, typename Number>
1383  const unsigned int geometry_index) const
1384  {
1385  return cell_type[geometry_index];
1386  }
1387 
1388 
1389 
1390  template <int dim, int spacedim, typename Number>
1393  const unsigned int cell_index) const
1394  {
1395  Assert(
1396  additional_data.store_cells,
1397  ExcMessage(
1398  "Cells have been not stored. You can enable this by Additional::store_cells."));
1399  return {triangulation.get(),
1400  cell_level_and_indices[cell_index].first,
1401  cell_level_and_indices[cell_index].second};
1402  }
1403 
1404 
1405 
1406  template <int dim, int spacedim, typename Number>
1407  template <typename NumberType>
1408  unsigned int
1410  const unsigned int n_q_points_unvectorized)
1411  {
1412  const unsigned int n_lanes =
1414  const unsigned int n_filled_lanes_last_batch =
1415  n_q_points_unvectorized % n_lanes;
1416  unsigned int n_q_points = n_q_points_unvectorized / n_lanes;
1417  if (n_filled_lanes_last_batch > 0)
1418  ++n_q_points;
1419  return n_q_points;
1420  }
1421 
1422 
1423 
1424  template <int dim, int spacedim, typename Number>
1425  template <bool is_face>
1426  unsigned int
1428  const unsigned int cell_index,
1429  const unsigned int face_number) const
1430  {
1432  return 0;
1433 
1434  const unsigned int compressed_cell_index =
1435  compute_compressed_cell_index(cell_index);
1436  if (!is_face)
1437  {
1438  Assert(state == State::cell_vector,
1439  ExcMessage(
1440  "This mapping info is not reinitialized for a cell vector!"));
1441  return compressed_cell_index;
1442  }
1443  else
1444  {
1446  ExcMessage(
1447  "cell_index has to be set if face_number is specified!"));
1448  Assert(state == State::faces_on_cells_in_vector,
1449  ExcMessage("This mapping info is not reinitialized for faces"
1450  " on cells in a vector!"));
1451  return cell_index_offset[compressed_cell_index] + face_number;
1452  }
1453  }
1454 
1455 
1456 
1457  template <int dim, int spacedim, typename Number>
1458  unsigned int
1460  const unsigned int cell_index) const
1461  {
1462  if (do_cell_index_compression)
1463  {
1464  Assert(cell_index_to_compressed_cell_index[cell_index] !=
1466  ExcMessage("Mapping info object was not initialized for this"
1467  " active cell index!"));
1468  return cell_index_to_compressed_cell_index[cell_index];
1469  }
1470  else
1471  return cell_index;
1472  }
1473 
1474 
1475  template <int dim, int spacedim, typename Number>
1476  void
1478  const unsigned int unit_points_index_offset,
1479  const unsigned int n_q_points,
1480  const unsigned int n_q_points_unvectorized,
1481  const std::vector<Point<dim>> &points)
1482  {
1483  const unsigned int n_lanes =
1485 
1486  for (unsigned int q = 0; q < n_q_points; ++q)
1487  {
1488  const unsigned int offset = unit_points_index_offset + q;
1489  for (unsigned int v = 0;
1490  v < n_lanes && q * n_lanes + v < n_q_points_unvectorized;
1491  ++v)
1492  for (unsigned int d = 0; d < dim; ++d)
1494  unit_points[offset][d], v) = points[q * n_lanes + v][d];
1495  }
1496  }
1497 
1498 
1499 
1500  template <int dim, int spacedim, typename Number>
1501  void
1503  const unsigned int unit_points_index_offset,
1504  const unsigned int n_q_points,
1505  const unsigned int n_q_points_unvectorized,
1506  const std::vector<Point<dim - 1>> &points)
1507  {
1508  const unsigned int n_lanes =
1510 
1511  for (unsigned int q = 0; q < n_q_points; ++q)
1512  {
1513  const unsigned int offset = unit_points_index_offset + q;
1514  for (unsigned int v = 0;
1515  v < n_lanes && q * n_lanes + v < n_q_points_unvectorized;
1516  ++v)
1517  for (unsigned int d = 0; d < dim - 1; ++d)
1519  unit_points_faces[offset][d], v) = points[q * n_lanes + v][d];
1520  }
1521  }
1522 
1523 
1524 
1525  template <int dim, int spacedim, typename Number>
1526  void
1528  const unsigned int unit_points_index_offset,
1529  const unsigned int n_q_points,
1530  const unsigned int n_q_points_unvectorized,
1531  const MappingInfo::MappingData &mapping_data,
1532  const std::vector<double> &weights,
1533  const unsigned int compressed_unit_point_index_offset,
1534  const bool affine_cell)
1535  {
1536  const unsigned int n_lanes =
1538 
1539  for (unsigned int q = 0; q < n_q_points; ++q)
1540  {
1541  const unsigned int offset = unit_points_index_offset + q;
1542  const unsigned int compressed_offset =
1543  compressed_unit_point_index_offset + q;
1544  for (unsigned int v = 0;
1545  v < n_lanes && q * n_lanes + v < n_q_points_unvectorized;
1546  ++v)
1547  {
1548  if (q == 0 || !affine_cell)
1549  {
1550  if (update_flags_mapping & UpdateFlags::update_jacobians)
1551  for (unsigned int d = 0; d < dim; ++d)
1552  for (unsigned int s = 0; s < spacedim; ++s)
1554  jacobians[compressed_offset][d][s], v) =
1555  mapping_data.jacobians[q * n_lanes + v][d][s];
1556  if (update_flags_mapping &
1558  for (unsigned int d = 0; d < dim; ++d)
1559  for (unsigned int s = 0; s < spacedim; ++s)
1561  inverse_jacobians[compressed_offset][s][d], v) =
1562  mapping_data.inverse_jacobians[q * n_lanes + v][s][d];
1563  }
1564  if (update_flags_mapping & UpdateFlags::update_JxW_values)
1565  {
1566  if (additional_data.use_global_weights)
1567  {
1569  JxW_values[offset], v) = weights[q * n_lanes + v];
1570  }
1571  else
1572  {
1574  JxW_values[offset], v) =
1575  mapping_data.JxW_values[q * n_lanes + v];
1576  }
1577  }
1578  if (update_flags_mapping & UpdateFlags::update_normal_vectors)
1579  for (unsigned int s = 0; s < spacedim; ++s)
1581  normal_vectors[offset][s], v) =
1582  mapping_data.normal_vectors[q * n_lanes + v][s];
1583  if (update_flags_mapping & UpdateFlags::update_quadrature_points)
1584  for (unsigned int s = 0; s < spacedim; ++s)
1586  real_points[offset][s], v) =
1587  mapping_data.quadrature_points[q * n_lanes + v][s];
1588  }
1589  }
1590  }
1591 
1592 
1593 
1594  template <int dim, int spacedim, typename Number>
1595  void
1597  const unsigned int n_unit_point_batches)
1598  {
1599  unit_points.resize(n_unit_point_batches);
1600  }
1601 
1602 
1603 
1604  template <int dim, int spacedim, typename Number>
1605  void
1607  const unsigned int n_unit_point_batches)
1608  {
1609  unit_points_faces.resize(n_unit_point_batches);
1610  }
1611 
1612 
1613 
1614  template <int dim, int spacedim, typename Number>
1615  void
1617  const unsigned int n_data_point_batches)
1618  {
1619  if (update_flags_mapping & UpdateFlags::update_jacobians)
1620  jacobians.resize(n_data_point_batches);
1621  if (update_flags_mapping & UpdateFlags::update_inverse_jacobians)
1622  inverse_jacobians.resize(n_data_point_batches);
1623  if (update_flags_mapping & UpdateFlags::update_JxW_values)
1624  JxW_values.resize(n_data_point_batches);
1625  if (update_flags_mapping & UpdateFlags::update_normal_vectors)
1626  normal_vectors.resize(n_data_point_batches);
1627  if (update_flags_mapping & UpdateFlags::update_quadrature_points)
1628  real_points.resize(n_data_point_batches);
1629  }
1630 
1631 
1632 
1633  template <int dim, int spacedim, typename Number>
1634  inline const Point<
1635  dim,
1638  const unsigned int offset) const
1639  {
1640  return unit_points.data() + offset;
1641  }
1642 
1643 
1644 
1645  template <int dim, int spacedim, typename Number>
1646  inline const Point<
1647  dim - 1,
1650  const unsigned int offset) const
1651  {
1652  return unit_points_faces.data() + offset;
1653  }
1654 
1655 
1656 
1657  template <int dim, int spacedim, typename Number>
1658  inline const Point<dim, Number> *
1660  const unsigned int offset) const
1661  {
1662  return real_points.data() + offset;
1663  }
1664 
1665 
1666 
1667  template <int dim, int spacedim, typename Number>
1668  unsigned int
1670  const unsigned int geometry_index) const
1671  {
1672  return unit_points_index[geometry_index];
1673  }
1674 
1675 
1676 
1677  template <int dim, int spacedim, typename Number>
1678  unsigned int
1680  const unsigned int geometry_index) const
1681  {
1682  return data_index_offsets[geometry_index];
1683  }
1684 
1685 
1686  template <int dim, int spacedim, typename Number>
1687  unsigned int
1689  const unsigned int geometry_index) const
1690  {
1691  return compressed_data_index_offsets[geometry_index];
1692  }
1693 
1694 
1695 
1696  template <int dim, int spacedim, typename Number>
1699  const unsigned int offset) const
1700  {
1701  return jacobians.data() + offset;
1702  }
1703 
1704 
1705 
1706  template <int dim, int spacedim, typename Number>
1709  const unsigned int offset) const
1710  {
1711  return inverse_jacobians.data() + offset;
1712  }
1713 
1714 
1715  template <int dim, int spacedim, typename Number>
1716  inline const Tensor<1, spacedim, Number> *
1718  const unsigned int offset) const
1719  {
1720  return normal_vectors.data() + offset;
1721  }
1722 
1723 
1724 
1725  template <int dim, int spacedim, typename Number>
1726  inline const Number *
1727  MappingInfo<dim, spacedim, Number>::get_JxW(const unsigned int offset) const
1728  {
1729  return JxW_values.data() + offset;
1730  }
1731 
1732 
1733 
1734  template <int dim, int spacedim, typename Number>
1735  const Mapping<dim, spacedim> &
1737  {
1738  return *mapping;
1739  }
1740 
1741 
1742 
1743  template <int dim, int spacedim, typename Number>
1744  UpdateFlags
1746  {
1747  return update_flags;
1748  }
1749 
1750 
1751 
1752  template <int dim, int spacedim, typename Number>
1753  UpdateFlags
1755  {
1756  return update_flags_mapping;
1757  }
1758 
1759 
1760 
1761  template <int dim, int spacedim, typename Number>
1762  boost::signals2::connection
1764  const std::function<void()> &set_is_reinitialized)
1765  {
1766  return is_reinitialized.connect(set_is_reinitialized);
1767  }
1768 
1769 
1770 
1771  template <int dim, int spacedim, typename Number>
1772  std::size_t
1774  {
1775  std::size_t memory = MemoryConsumption::memory_consumption(unit_points);
1776  memory += MemoryConsumption::memory_consumption(unit_points_faces);
1777  memory += MemoryConsumption::memory_consumption(unit_points_index);
1778  memory += cell_type.capacity() *
1780  memory += MemoryConsumption::memory_consumption(data_index_offsets);
1781  memory +=
1782  MemoryConsumption::memory_consumption(compressed_data_index_offsets);
1783  memory += MemoryConsumption::memory_consumption(JxW_values);
1784  memory += MemoryConsumption::memory_consumption(normal_vectors);
1785  memory += MemoryConsumption::memory_consumption(jacobians);
1786  memory += MemoryConsumption::memory_consumption(inverse_jacobians);
1787  memory += MemoryConsumption::memory_consumption(real_points);
1788  memory += MemoryConsumption::memory_consumption(n_q_points_unvectorized);
1789  memory += MemoryConsumption::memory_consumption(cell_index_offset);
1791  cell_index_to_compressed_cell_index);
1792  memory += MemoryConsumption::memory_consumption(cell_level_and_indices);
1793  memory += sizeof(*this);
1794  return memory;
1795  }
1796 } // namespace NonMatching
1797 
1799 
1800 #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
const UpdateFlags update_flags
Definition: mapping_info.h:650
unsigned int compute_compressed_data_index_offset(const unsigned int geometry_index) const
AlignedVector< Point< dim, VectorizedArrayType > > unit_points
Definition: mapping_info.h:622
const DerivativeForm< 1, spacedim, dim, Number > * get_inverse_jacobian(const unsigned int offset) const
AlignedVector< Number > JxW_values
Definition: mapping_info.h:686
::internal::MatrixFreeFunctions::GeometryType get_cell_type(const unsigned int geometry_index) const
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)
std::shared_ptr< typename Mapping< dim, spacedim >::InternalDataBase > internal_mapping_data
Definition: mapping_info.h:640
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
Definition: mapping_info.h:701
SmartPointer< const Triangulation< dim, spacedim > > triangulation
Definition: mapping_info.h:752
const Point< dim, VectorizedArrayType > * get_unit_point(const unsigned int offset) const
unsigned int get_n_q_points_unvectorized(const unsigned int geometry_index) 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
Definition: mapping_info.h:728
MappingInfo(const Mapping< dim > &mapping, const UpdateFlags update_flags, const AdditionalData additional_data=AdditionalData())
Definition: mapping_info.h:765
const SmartPointer< const Mapping< dim, spacedim > > mapping
Definition: mapping_info.h:645
AlignedVector< Tensor< 1, spacedim, Number > > normal_vectors
Definition: mapping_info.h:693
unsigned int compute_geometry_index_offset(const unsigned int cell_index, const unsigned int face_number) const
std::vector<::internal::MatrixFreeFunctions::GeometryType > cell_type
Definition: mapping_info.h:668
void reinit(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const std::vector< Point< dim >> &unit_points)
Definition: mapping_info.h:820
const DerivativeForm< 1, dim, spacedim, Number > * get_jacobian(const unsigned int offset) const
boost::signals2::signal< void()> is_reinitialized
Definition: mapping_info.h:745
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)
void resize_unit_points_faces(const unsigned int n_unit_point_batches)
const Tensor< 1, spacedim, Number > * get_normal_vector(const unsigned int offset) const
typename ::internal::VectorizedArrayTrait< Number >::vectorized_value_type VectorizedArrayType
Definition: mapping_info.h:254
void resize_unit_points(const unsigned int n_unit_point_batches)
std::size_t memory_consumption() const
UpdateFlags get_update_flags() const
void resize_data_fields(const unsigned int n_data_point_batches)
unsigned int compute_unit_point_index_offset(const unsigned int geometry_index) const
std::vector< unsigned int > cell_index_to_compressed_cell_index
Definition: mapping_info.h:734
AlignedVector< Point< spacedim, Number > > real_points
Definition: mapping_info.h:717
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, const unsigned int compressed_unit_point_index_offset, const bool affine_cell)
unsigned int compute_data_index_offset(const unsigned int geometry_index) const
AlignedVector< DerivativeForm< 1, spacedim, dim, Number > > inverse_jacobians
Definition: mapping_info.h:710
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)
Definition: mapping_info.h:915
std::vector< unsigned int > unit_points_index
Definition: mapping_info.h:634
std::vector< std::pair< int, int > > cell_level_and_indices
Definition: mapping_info.h:758
std::vector< unsigned int > compressed_data_index_offsets
Definition: mapping_info.h:678
std::vector< unsigned int > data_index_offsets
Definition: mapping_info.h:673
std::vector< unsigned int > n_q_points_unvectorized
Definition: mapping_info.h:722
const Number * get_JxW(const unsigned int offset) const
const AdditionalData additional_data
Definition: mapping_info.h:660
AlignedVector< Point< dim - 1, VectorizedArrayType > > unit_points_faces
Definition: mapping_info.h:629
unsigned int compute_n_q_points(const unsigned int n_q_points_unvectorized)
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 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)
UpdateFlags update_flags_mapping
Definition: mapping_info.h:655
void do_reinit_cells(const ContainerType &cell_iterator_range, const std::vector< QuadratureType > &quadrature_vector, const unsigned int n_unfiltered_cells, const std::function< void(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const QuadratureType &quadrature, MappingData &mapping_data)> &compute_mapping_data)
Definition: mapping_info.h:938
UpdateFlags get_update_flags_mapping() const
Triangulation< dim, spacedim >::cell_iterator get_cell_iterator(const unsigned int cell_index) const
static UpdateFlags required_update_flags(const SmartPointer< const Mapping< dim, spacedim >> mapping, const UpdateFlags &update_flags)
Definition: mapping_info.h:55
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)
Definition: mapping_info.h:139
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)
Definition: mapping_info.h:63
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)
Definition: mapping_info.h:103
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< Point< dim > > & get_points() const
const std::vector< double > & get_weights() const
bool empty() const
unsigned int size() const
Triangulation< dim, spacedim > & get_triangulation()
std::vector< DerivativeForm< 1, spacedim, dim > > inverse_jacobians
void initialize(const unsigned int n_quadrature_points, const UpdateFlags flags)
std::vector< Tensor< 1, spacedim > > normal_vectors
std::vector< DerivativeForm< 1, dim, spacedim > > jacobians
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
unsigned int cell_index
Definition: grid_tools.cc:1192
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
Definition: exceptions.h:1616
static ::ExceptionBase & ExcNotImplemented()
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1789
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.
Expression fabs(const Expression &x)
double diameter(const Triangulation< dim, spacedim > &tria)
Definition: grid_tools.cc:111
@ general
No special properties.
std::enable_if_t< std::is_fundamental_v< T >, std::size_t > memory_consumption(const T &t)
::internal::MatrixFreeFunctions::GeometryType compute_geometry_type(const double diameter, const std::vector< DerivativeForm< 1, dim, spacedim, double >> &inverse_jacobians)
Definition: mapping_info.h:188
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > epsilon(const Tensor< 2, dim, Number > &Grad_u)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
unsigned int n_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
Definition: tria.cc:14833
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618
static const unsigned int invalid_unsigned_int
Definition: types.h:213
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
ComparisonResult compare(const std::vector< T > &v1, const std::vector< T > &v2) const
AdditionalData(const bool use_global_weights=false, const bool store_cells=false)
Definition: mapping_info.h:265
static constexpr std::size_t width()
static value_type & get(value_type &value, unsigned int c)