Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
mapping_info.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2011 - 2019 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_matrix_free_mapping_info_h
18 #define dealii_matrix_free_mapping_info_h
19 
20 
21 #include <deal.II/base/aligned_vector.h>
22 #include <deal.II/base/exceptions.h>
23 #include <deal.II/base/vectorization.h>
24 
25 #include <deal.II/fe/fe.h>
26 #include <deal.II/fe/mapping.h>
27 
28 #include <deal.II/hp/q_collection.h>
29 
30 #include <deal.II/matrix_free/face_info.h>
31 #include <deal.II/matrix_free/helper_functions.h>
32 
33 #include <memory>
34 
35 
36 DEAL_II_NAMESPACE_OPEN
37 
38 
39 namespace internal
40 {
41  namespace MatrixFreeFunctions
42  {
50  enum GeometryType : unsigned char
51  {
55  cartesian = 0,
56 
60  affine = 1,
61 
67 
72  general = 3
73  };
74 
75 
76 
111  template <int structdim, int spacedim, typename Number>
113  {
114  struct QuadratureDescriptor
115  {
119  QuadratureDescriptor();
120 
124  void
125  initialize(const Quadrature<1> &quadrature_1d,
126  const UpdateFlags update_flags_inner_faces = update_default);
127 
131  std::size_t
132  memory_consumption() const;
133 
137  unsigned int n_q_points;
138 
142  Quadrature<structdim> quadrature;
143 
148  std::array<AlignedVector<Number>, structdim> tensor_quadrature_weights;
149 
155  AlignedVector<Number> quadrature_weights;
156 
163  ::Table<2, unsigned int> face_orientations;
164  };
165 
173  std::vector<QuadratureDescriptor> descriptor;
174 
182 
191 
199 
212 
229  spacedim *(spacedim + 1) / 2,
230  Tensor<1, spacedim, VectorizedArray<Number>>>>
232 
242 
250 
259 
266  unsigned int
267  quad_index_from_n_q_points(const unsigned int n_q_points) const;
268 
273  template <typename StreamType>
274  void
275  print_memory_consumption(StreamType & out,
276  const SizeInfo &task_info) const;
277 
281  std::size_t
282  memory_consumption() const;
283  };
284 
285 
286 
295  template <int dim, typename Number>
296  struct MappingInfo
297  {
306  void
307  initialize(
308  const ::Triangulation<dim> & tria,
309  const std::vector<std::pair<unsigned int, unsigned int>> & cells,
311  const std::vector<unsigned int> & active_fe_index,
312  const Mapping<dim> & mapping,
313  const std::vector<::hp::QCollection<1>> &quad,
314  const UpdateFlags update_flags_cells,
315  const UpdateFlags update_flags_boundary_faces,
316  const UpdateFlags update_flags_inner_faces,
317  const UpdateFlags update_flags_faces_by_cells);
318 
323  get_cell_type(const unsigned int cell_chunk_no) const;
324 
328  void
329  clear();
330 
334  std::size_t
335  memory_consumption() const;
336 
341  template <typename StreamType>
342  void
343  print_memory_consumption(StreamType & out,
344  const TaskInfo &task_info) const;
345 
352  std::vector<GeometryType> cell_type;
353 
361  std::vector<GeometryType> face_type;
362 
366  std::vector<MappingInfoStorage<dim, dim, Number>> cell_data;
367 
371  std::vector<MappingInfoStorage<dim - 1, dim, Number>> face_data;
372 
377  std::vector<MappingInfoStorage<dim - 1, dim, Number>> face_data_by_cells;
378 
383  void
385  const ::Triangulation<dim> & tria,
386  const std::vector<std::pair<unsigned int, unsigned int>> &cells,
387  const std::vector<unsigned int> & active_fe_index,
388  const Mapping<dim> & mapping,
389  const std::vector<::hp::QCollection<1>> &quad,
390  const UpdateFlags update_flags_cells);
391 
396  void
398  const ::Triangulation<dim> & tria,
399  const std::vector<std::pair<unsigned int, unsigned int>> &cells,
400  const std::vector<
402  const Mapping<dim> & mapping,
403  const std::vector<::hp::QCollection<1>> &quad,
404  const UpdateFlags update_flags_boundary_faces,
405  const UpdateFlags update_flags_inner_faces);
406 
411  void
413  const ::Triangulation<dim> & tria,
414  const std::vector<std::pair<unsigned int, unsigned int>> &cells,
415  const Mapping<dim> & mapping,
416  const std::vector<::hp::QCollection<1>> & quad,
417  const UpdateFlags update_flags_faces_by_cells);
418 
423  static UpdateFlags
424  compute_update_flags(const UpdateFlags update_flags,
425  const std::vector<::hp::QCollection<1>> &quad =
426  std::vector<::hp::QCollection<1>>());
427  };
428 
429 
430 
437  template <int, typename, bool>
439 
440  template <int dim, typename Number>
441  struct MappingInfoCellsOrFaces<dim, Number, false>
442  {
444  get(const MappingInfo<dim, Number> &mapping_info,
445  const unsigned int quad_no)
446  {
447  AssertIndexRange(quad_no, mapping_info.cell_data.size());
448  return &mapping_info.cell_data[quad_no];
449  }
450  };
451 
452  template <int dim, typename Number>
453  struct MappingInfoCellsOrFaces<dim, Number, true>
454  {
455  static const MappingInfoStorage<dim - 1, dim, Number> *
456  get(const MappingInfo<dim, Number> &mapping_info,
457  const unsigned int quad_no)
458  {
459  AssertIndexRange(quad_no, mapping_info.face_data.size());
460  return &mapping_info.face_data[quad_no];
461  }
462  };
463 
464 
465 
477  template <typename Number>
479  {
480  FPArrayComparator(const Number scaling);
481 
482  bool
483  operator()(const std::vector<Number> &v1,
484  const std::vector<Number> &v2) const;
485 
486  bool
487  operator()(
490  const;
491 
492  template <int dim>
493  bool
494  operator()(
495  const Tensor<
496  1,
497  dim,
499  const Tensor<
500  1,
501  dim,
503  const;
504 
505  template <int dim>
506  bool
507  operator()(
508  const Tensor<
509  2,
510  dim,
512  const Tensor<
513  2,
514  dim,
516  const;
517 
518  Number tolerance;
519  };
520 
521 
522 
523  /* ------------------- inline functions ----------------------------- */
524 
525  template <int structdim, int spacedim, typename Number>
526  inline unsigned int
528  const unsigned int n_q_points) const
529  {
530  for (unsigned int i = 0; i < descriptor.size(); ++i)
531  if (n_q_points == descriptor[i].n_q_points)
532  return i;
533  return 0;
534  }
535 
536 
537 
538  template <int dim, typename Number>
539  inline GeometryType
540  MappingInfo<dim, Number>::get_cell_type(const unsigned int cell_no) const
541  {
542  AssertIndexRange(cell_no, cell_type.size());
543  return cell_type[cell_no];
544  }
545 
546  } // end of namespace MatrixFreeFunctions
547 } // end of namespace internal
548 
549 DEAL_II_NAMESPACE_CLOSE
550 
551 #endif
AlignedVector< Tensor< 1, spacedim, VectorizedArray< Number > > > normals_times_jacobians[2]
Definition: mapping_info.h:241
AlignedVector< unsigned int > data_index_offsets
Definition: mapping_info.h:181
AlignedVector< Tensor< 1, spacedim *(spacedim+1)/2, Tensor< 1, spacedim, VectorizedArray< Number > > > > jacobian_gradients[2]
Definition: mapping_info.h:231
AlignedVector< unsigned int > quadrature_point_offsets
Definition: mapping_info.h:249
#define AssertIndexRange(index, range)
Definition: exceptions.h:1637
void initialize(const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int >> &cells, const FaceInfo< VectorizedArray< Number >::n_array_elements > &faces, const std::vector< unsigned int > &active_fe_index, const Mapping< dim > &mapping, const std::vector<::hp::QCollection< 1 >> &quad, const UpdateFlags update_flags_cells, const UpdateFlags update_flags_boundary_faces, const UpdateFlags update_flags_inner_faces, const UpdateFlags update_flags_faces_by_cells)
void initialize_cells(const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int >> &cells, const std::vector< unsigned int > &active_fe_index, const Mapping< dim > &mapping, const std::vector<::hp::QCollection< 1 >> &quad, const UpdateFlags update_flags_cells)
void print_memory_consumption(StreamType &out, const TaskInfo &task_info) const
std::vector< MappingInfoStorage< dim, dim, Number > > cell_data
Definition: mapping_info.h:366
std::vector< MappingInfoStorage< dim - 1, dim, Number > > face_data_by_cells
Definition: mapping_info.h:377
std::vector< GeometryType > face_type
Definition: mapping_info.h:361
std::vector< QuadratureDescriptor > descriptor
Definition: mapping_info.h:173
No update.
UpdateFlags
GeometryType get_cell_type(const unsigned int cell_chunk_no) const
Definition: mapping_info.h:540
void initialize_faces_by_cells(const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int >> &cells, const Mapping< dim > &mapping, const std::vector<::hp::QCollection< 1 >> &quad, const UpdateFlags update_flags_faces_by_cells)
std::vector< MappingInfoStorage< dim - 1, dim, Number > > face_data
Definition: mapping_info.h:371
void print_memory_consumption(StreamType &out, const SizeInfo &task_info) const
void initialize_faces(const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int >> &cells, const std::vector< FaceToCellTopology< VectorizedArray< Number >::n_array_elements >> &faces, const Mapping< dim > &mapping, const std::vector<::hp::QCollection< 1 >> &quad, const UpdateFlags update_flags_boundary_faces, const UpdateFlags update_flags_inner_faces)
AlignedVector< Point< spacedim, VectorizedArray< Number > > > quadrature_points
Definition: mapping_info.h:258
unsigned int quad_index_from_n_q_points(const unsigned int n_q_points) const
Definition: mapping_info.h:527
Definition: mpi.h:90
AlignedVector< Tensor< 1, spacedim, VectorizedArray< Number > > > normal_vectors
Definition: mapping_info.h:198
static UpdateFlags compute_update_flags(const UpdateFlags update_flags, const std::vector<::hp::QCollection< 1 >> &quad=std::vector<::hp::QCollection< 1 >>())
std::vector< GeometryType > cell_type
Definition: mapping_info.h:352
AlignedVector< Tensor< 2, spacedim, VectorizedArray< Number > > > jacobians[2]
Definition: mapping_info.h:211
AlignedVector< VectorizedArray< Number > > JxW_values
Definition: mapping_info.h:190