Reference documentation for deal.II version 9.3.3
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
mapping_info.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2011 - 2021 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16
17#ifndef dealii_matrix_free_mapping_info_h
18#define dealii_matrix_free_mapping_info_h
19
20
21#include <deal.II/base/config.h>
22
26
27#include <deal.II/fe/fe.h>
28#include <deal.II/fe/mapping.h>
29
31
34
37
38#include <memory>
39
40
42
43
44namespace internal
45{
46 namespace MatrixFreeFunctions
47 {
55 enum GeometryType : unsigned char
56 {
61
65 affine = 1,
66
72
77 general = 3
78 };
79
80
81
114 template <int structdim,
115 int spacedim,
116 typename Number,
117 typename VectorizedArrayType>
119 {
120 static_assert(
121 std::is_same<Number, typename VectorizedArrayType::value_type>::value,
122 "Type of Number and of VectorizedArrayType do not match.");
123
125 {
130
134 template <int dim_q>
135 void
137 const UpdateFlags update_flags_inner_faces = update_default);
138
142 void
144 const UpdateFlags update_flags_inner_faces = update_default);
145
149 std::size_t
151
155 unsigned int n_q_points;
156
162
167
172 std::array<AlignedVector<Number>, structdim> tensor_quadrature_weights;
173
180
188 };
189
197 std::vector<QuadratureDescriptor> descriptor;
198
205 std::vector<::hp::QCollection<structdim>> q_collection;
206
214
223
230
242 std::array<AlignedVector<Tensor<2, spacedim, VectorizedArrayType>>, 2>
244
260 std::array<
262 spacedim *(spacedim + 1) / 2,
264 2>
266
274 std::array<AlignedVector<Tensor<1, spacedim, VectorizedArrayType>>, 2>
276
284
293
297 void
299
306 unsigned int
307 quad_index_from_n_q_points(const unsigned int n_q_points) const;
308
313 template <typename StreamType>
314 void
315 print_memory_consumption(StreamType & out,
316 const TaskInfo &task_info) const;
317
321 std::size_t
323 };
324
325
326
333 template <int dim, typename Number, typename VectorizedArrayType>
335 {
344 void
346 const ::Triangulation<dim> & tria,
347 const std::vector<std::pair<unsigned int, unsigned int>> &cells,
348 const FaceInfo<VectorizedArrayType::size()> & faces,
349 const std::vector<unsigned int> &active_fe_index,
350 const std::shared_ptr<::hp::MappingCollection<dim>> &mapping,
351 const std::vector<::hp::QCollection<dim>> & quad,
356
363 void
365 const ::Triangulation<dim> & tria,
366 const std::vector<std::pair<unsigned int, unsigned int>> &cells,
367 const FaceInfo<VectorizedArrayType::size()> & faces,
368 const std::vector<unsigned int> &active_fe_index,
369 const std::shared_ptr<::hp::MappingCollection<dim>> &mapping);
370
375 get_cell_type(const unsigned int cell_chunk_no) const;
376
380 void
382
386 std::size_t
388
393 template <typename StreamType>
394 void
395 print_memory_consumption(StreamType & out,
396 const TaskInfo &task_info) const;
397
402
408
414
420
427 std::vector<GeometryType> cell_type;
428
436 std::vector<GeometryType> face_type;
437
441 std::vector<MappingInfoStorage<dim, dim, Number, VectorizedArrayType>>
443
447 std::vector<MappingInfoStorage<dim - 1, dim, Number, VectorizedArrayType>>
449
454 std::vector<MappingInfoStorage<dim - 1, dim, Number, VectorizedArrayType>>
456
460 std::shared_ptr<::hp::MappingCollection<dim>> mapping_collection;
461
466
471 std::vector<std::vector<::ReferenceCell>> reference_cell_types;
472
494 void
496 const ::Triangulation<dim> & tria,
497 const std::vector<std::pair<unsigned int, unsigned int>> &cells,
498 const std::vector<FaceToCellTopology<VectorizedArrayType::size()>>
499 &faces);
500
505 void
507 const ::Triangulation<dim> & tria,
508 const std::vector<std::pair<unsigned int, unsigned int>> &cells,
509 const std::vector<unsigned int> & active_fe_index,
510 const ::hp::MappingCollection<dim> &mapping);
511
516 void
518 const ::Triangulation<dim> & tria,
519 const std::vector<std::pair<unsigned int, unsigned int>> &cells,
520 const std::vector<FaceToCellTopology<VectorizedArrayType::size()>>
521 & faces,
522 const std::vector<unsigned int> & active_fe_index,
523 const ::hp::MappingCollection<dim> &mapping);
524
529 void
531 const ::Triangulation<dim> & tria,
532 const std::vector<std::pair<unsigned int, unsigned int>> &cells,
533 const ::hp::MappingCollection<dim> & mapping);
534
539 static UpdateFlags
541 const UpdateFlags update_flags,
542 const std::vector<::hp::QCollection<dim>> &quad =
543 std::vector<::hp::QCollection<dim>>());
544 };
545
546
547
552 template <int, typename, bool, typename>
554
555 template <int dim, typename Number, typename VectorizedArrayType>
556 struct MappingInfoCellsOrFaces<dim, Number, false, VectorizedArrayType>
557 {
560 const unsigned int quad_no)
561 {
562 AssertIndexRange(quad_no, mapping_info.cell_data.size());
563 return &mapping_info.cell_data[quad_no];
564 }
565 };
566
567 template <int dim, typename Number, typename VectorizedArrayType>
568 struct MappingInfoCellsOrFaces<dim, Number, true, VectorizedArrayType>
569 {
570 static const MappingInfoStorage<dim - 1, dim, Number, VectorizedArrayType>
571 *
573 const unsigned int quad_no)
574 {
575 AssertIndexRange(quad_no, mapping_info.face_data.size());
576 return &mapping_info.face_data[quad_no];
577 }
578 };
579
580
581
593 template <typename Number,
594 typename VectorizedArrayType = VectorizedArray<Number>>
596 {
597 FPArrayComparator(const Number scaling);
598
602 bool
603 operator()(const std::vector<Number> &v1,
604 const std::vector<Number> &v2) const;
605
610 bool
612 const Tensor<1, VectorizedArrayType::size(), Number> &t1,
613 const Tensor<1, VectorizedArrayType::size(), Number> &t2) const;
614
619 template <int dim>
620 bool
622 const Tensor<1, dim, Tensor<1, VectorizedArrayType::size(), Number>>
623 &t1,
624 const Tensor<1, dim, Tensor<1, VectorizedArrayType::size(), Number>>
625 &t2) const;
626
631 template <int dim>
632 bool
634 const Tensor<2, dim, Tensor<1, VectorizedArrayType::size(), Number>>
635 &t1,
636 const Tensor<2, dim, Tensor<1, VectorizedArrayType::size(), Number>>
637 &t2) const;
638
642 template <int dim>
643 bool
644 operator()(const std::array<Tensor<2, dim, Number>, dim + 1> &t1,
645 const std::array<Tensor<2, dim, Number>, dim + 1> &t2) const;
646
647 Number tolerance;
648 };
649
650
651
652 /* ------------------- inline functions ----------------------------- */
653
654 template <int structdim,
655 int spacedim,
656 typename Number,
657 typename VectorizedArrayType>
658 inline unsigned int
660 quad_index_from_n_q_points(const unsigned int n_q_points) const
661 {
662 for (unsigned int i = 0; i < descriptor.size(); ++i)
663 if (n_q_points == descriptor[i].n_q_points)
664 return i;
665 return 0;
666 }
667
668
669
670 template <int dim, typename Number, typename VectorizedArrayType>
671 inline GeometryType
673 const unsigned int cell_no) const
674 {
675 AssertIndexRange(cell_no, cell_type.size());
676 return cell_type[cell_no];
677 }
678
679 } // end of namespace MatrixFreeFunctions
680} // end of namespace internal
681
683
684#endif
Definition: tensor.h:472
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:402
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:403
UpdateFlags
@ update_default
No update.
const unsigned int v1
Definition: grid_tools.cc:963
#define AssertIndexRange(index, range)
Definition: exceptions.h:1690
bool operator()(const std::vector< Number > &v1, const std::vector< Number > &v2) const
bool operator()(const std::array< Tensor< 2, dim, Number >, dim+1 > &t1, const std::array< Tensor< 2, dim, Number >, dim+1 > &t2) const
bool operator()(const Tensor< 1, dim, Tensor< 1, VectorizedArrayType::size(), Number > > &t1, const Tensor< 1, dim, Tensor< 1, VectorizedArrayType::size(), Number > > &t2) const
bool operator()(const Tensor< 1, VectorizedArrayType::size(), Number > &t1, const Tensor< 1, VectorizedArrayType::size(), Number > &t2) const
bool operator()(const Tensor< 2, dim, Tensor< 1, VectorizedArrayType::size(), Number > > &t1, const Tensor< 2, dim, Tensor< 1, VectorizedArrayType::size(), Number > > &t2) const
static const MappingInfoStorage< dim, dim, Number, VectorizedArrayType > * get(const MappingInfo< dim, Number, VectorizedArrayType > &mapping_info, const unsigned int quad_no)
Definition: mapping_info.h:559
static const MappingInfoStorage< dim - 1, dim, Number, VectorizedArrayType > * get(const MappingInfo< dim, Number, VectorizedArrayType > &mapping_info, const unsigned int quad_no)
Definition: mapping_info.h:572
void initialize(const Quadrature< dim_q > &quadrature, const UpdateFlags update_flags_inner_faces=update_default)
void initialize(const Quadrature< 1 > &quadrature_1d, const UpdateFlags update_flags_inner_faces=update_default)
std::array< AlignedVector< Number >, structdim > tensor_quadrature_weights
Definition: mapping_info.h:172
std::array< AlignedVector< Tensor< 1, spacedim *(spacedim+1)/2, Tensor< 1, spacedim, VectorizedArrayType > > >, 2 > jacobian_gradients
Definition: mapping_info.h:265
AlignedVector< Point< spacedim, VectorizedArrayType > > quadrature_points
Definition: mapping_info.h:292
unsigned int quad_index_from_n_q_points(const unsigned int n_q_points) const
Definition: mapping_info.h:660
std::vector<::hp::QCollection< structdim > > q_collection
Definition: mapping_info.h:205
AlignedVector< unsigned int > quadrature_point_offsets
Definition: mapping_info.h:283
AlignedVector< Tensor< 1, spacedim, VectorizedArrayType > > normal_vectors
Definition: mapping_info.h:229
AlignedVector< VectorizedArrayType > JxW_values
Definition: mapping_info.h:222
AlignedVector< unsigned int > data_index_offsets
Definition: mapping_info.h:213
std::array< AlignedVector< Tensor< 1, spacedim, VectorizedArrayType > >, 2 > normals_times_jacobians
Definition: mapping_info.h:275
void print_memory_consumption(StreamType &out, const TaskInfo &task_info) const
std::array< AlignedVector< Tensor< 2, spacedim, VectorizedArrayType > >, 2 > jacobians
Definition: mapping_info.h:243
std::vector< QuadratureDescriptor > descriptor
Definition: mapping_info.h:197
std::vector< MappingInfoStorage< dim, dim, Number, VectorizedArrayType > > cell_data
Definition: mapping_info.h:442
void initialize_faces(const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int > > &cells, const std::vector< FaceToCellTopology< VectorizedArrayType::size()> > &faces, const std::vector< unsigned int > &active_fe_index, const ::hp::MappingCollection< dim > &mapping)
std::vector< MappingInfoStorage< dim - 1, dim, Number, VectorizedArrayType > > face_data_by_cells
Definition: mapping_info.h:455
std::shared_ptr<::hp::MappingCollection< dim > > mapping_collection
Definition: mapping_info.h:460
GeometryType get_cell_type(const unsigned int cell_chunk_no) const
Definition: mapping_info.h:672
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 ::hp::MappingCollection< dim > &mapping)
std::vector< GeometryType > cell_type
Definition: mapping_info.h:427
void initialize_faces_by_cells(const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int > > &cells, const ::hp::MappingCollection< dim > &mapping)
std::vector< GeometryType > face_type
Definition: mapping_info.h:436
std::vector< MappingInfoStorage< dim - 1, dim, Number, VectorizedArrayType > > face_data
Definition: mapping_info.h:448
std::vector< std::vector<::ReferenceCell > > reference_cell_types
Definition: mapping_info.h:471
SmartPointer< const Mapping< dim > > mapping
Definition: mapping_info.h:465
void initialize(const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int > > &cells, const FaceInfo< VectorizedArrayType::size()> &faces, const std::vector< unsigned int > &active_fe_index, const std::shared_ptr<::hp::MappingCollection< dim > > &mapping, const std::vector<::hp::QCollection< dim > > &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 update_mapping(const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int > > &cells, const FaceInfo< VectorizedArrayType::size()> &faces, const std::vector< unsigned int > &active_fe_index, const std::shared_ptr<::hp::MappingCollection< dim > > &mapping)
void print_memory_consumption(StreamType &out, const TaskInfo &task_info) const
static UpdateFlags compute_update_flags(const UpdateFlags update_flags, const std::vector<::hp::QCollection< dim > > &quad=std::vector<::hp::QCollection< dim > >())
void compute_mapping_q(const ::Triangulation< dim > &tria, const std::vector< std::pair< unsigned int, unsigned int > > &cells, const std::vector< FaceToCellTopology< VectorizedArrayType::size()> > &faces)