Reference documentation for deal.II version GIT 49863513f1 2023-10-02 19:40: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_q.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2000 - 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 #ifndef dealii_mapping_q_h
17 #define dealii_mapping_q_h
18 
19 
20 #include <deal.II/base/config.h>
21 
25 #include <deal.II/base/table.h>
27 
28 #include <deal.II/fe/mapping.h>
29 
31 
33 
34 #include <array>
35 #include <cmath>
36 
38 
39 #ifndef DOXYGEN
40 template <int, int>
41 class MappingQCache;
42 #endif
43 
109 template <int dim, int spacedim = dim>
110 class MappingQ : public Mapping<dim, spacedim>
111 {
112 public:
118  MappingQ(const unsigned int polynomial_degree);
119 
126  MappingQ(const unsigned int polynomial_degree,
127  const bool use_mapping_q_on_all_cells);
128 
132  MappingQ(const MappingQ<dim, spacedim> &mapping);
133 
134  // for documentation, see the Mapping base class
135  virtual std::unique_ptr<Mapping<dim, spacedim>>
136  clone() const override;
137 
142  unsigned int
143  get_degree() const;
144 
149  virtual bool
150  preserves_vertex_locations() const override;
151 
152  // for documentation, see the Mapping base class
153  virtual BoundingBox<spacedim>
155  &cell) const override;
156 
157  virtual bool
158  is_compatible_with(const ReferenceCell &reference_cell) const override;
159 
165  // for documentation, see the Mapping base class
166  virtual Point<spacedim>
168  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
169  const Point<dim> &p) const override;
170 
171  // for documentation, see the Mapping base class
172  virtual Point<dim>
174  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
175  const Point<spacedim> &p) const override;
176 
177  // for documentation, see the Mapping base class
178  virtual void
180  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
181  const ArrayView<const Point<spacedim>> &real_points,
182  const ArrayView<Point<dim>> &unit_points) const override;
183 
193  // for documentation, see the Mapping base class
194  virtual void
195  transform(const ArrayView<const Tensor<1, dim>> &input,
196  const MappingKind kind,
197  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
198  const ArrayView<Tensor<1, spacedim>> &output) const override;
199 
200  // for documentation, see the Mapping base class
201  virtual void
203  const MappingKind kind,
204  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
205  const ArrayView<Tensor<2, spacedim>> &output) const override;
206 
207  // for documentation, see the Mapping base class
208  virtual void
209  transform(const ArrayView<const Tensor<2, dim>> &input,
210  const MappingKind kind,
211  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
212  const ArrayView<Tensor<2, spacedim>> &output) const override;
213 
214  // for documentation, see the Mapping base class
215  virtual void
217  const MappingKind kind,
218  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
219  const ArrayView<Tensor<3, spacedim>> &output) const override;
220 
221  // for documentation, see the Mapping base class
222  virtual void
223  transform(const ArrayView<const Tensor<3, dim>> &input,
224  const MappingKind kind,
225  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
226  const ArrayView<Tensor<3, spacedim>> &output) const override;
227 
249  void
251  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
252  const ArrayView<const Point<dim>> &unit_points,
253  const UpdateFlags update_flags,
255  &output_data) const;
256 
279  void
281  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
282  const unsigned int face_number,
283  const Quadrature<dim - 1> &face_quadrature,
284  const typename Mapping<dim, spacedim>::InternalDataBase &internal_data,
286  &output_data) const;
287 
304  class InternalData : public Mapping<dim, spacedim>::InternalDataBase
305  {
306  public:
311  InternalData(const unsigned int polynomial_degree);
312 
321  void
322  initialize(const UpdateFlags update_flags,
323  const Quadrature<dim> &quadrature,
324  const unsigned int n_original_q_points);
325 
331  void
332  initialize_face(const UpdateFlags update_flags,
333  const Quadrature<dim> &quadrature,
334  const unsigned int n_original_q_points);
335 
339  virtual std::size_t
340  memory_consumption() const override;
341 
347  std::vector<Point<dim>> quadrature_points;
348 
362  std::array<std::vector<Tensor<1, dim>>,
365 
370  const unsigned int polynomial_degree;
371 
381  const unsigned int n_shape_functions;
382 
383  /*
384  * The default line support points. Is used in when the shape function
385  * values are computed.
386  *
387  * The number of quadrature points depends on the degree of this
388  * class, and it matches the number of degrees of freedom of an
389  * FE_Q<1>(this->degree).
390  */
392 
409  std::min<std::size_t>(VectorizedArray<double>::size(),
410  (dim <= 2 ? 2 : 4))>;
411 
418 
424 
429 
433  mutable std::vector<AlignedVector<Tensor<1, spacedim>>> aux;
434 
439  mutable std::vector<Point<spacedim>> mapping_support_points;
440 
446 
452 
460  };
461 
462 protected:
463  // documentation can be found in Mapping::requires_update_flags()
464  virtual UpdateFlags
465  requires_update_flags(const UpdateFlags update_flags) const override;
466 
467  // documentation can be found in Mapping::get_data()
468  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
469  get_data(const UpdateFlags, const Quadrature<dim> &quadrature) const override;
470 
472 
473  // documentation can be found in Mapping::get_face_data()
474  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
475  get_face_data(const UpdateFlags flags,
476  const hp::QCollection<dim - 1> &quadrature) const override;
477 
478  // documentation can be found in Mapping::get_subface_data()
479  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
480  get_subface_data(const UpdateFlags flags,
481  const Quadrature<dim - 1> &quadrature) const override;
482 
483  // documentation can be found in Mapping::fill_fe_values()
486  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
487  const CellSimilarity::Similarity cell_similarity,
488  const Quadrature<dim> &quadrature,
489  const typename Mapping<dim, spacedim>::InternalDataBase &internal_data,
491  &output_data) const override;
492 
494 
495  // documentation can be found in Mapping::fill_fe_face_values()
496  virtual void
498  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
499  const unsigned int face_no,
500  const hp::QCollection<dim - 1> &quadrature,
501  const typename Mapping<dim, spacedim>::InternalDataBase &internal_data,
503  &output_data) const override;
504 
505  // documentation can be found in Mapping::fill_fe_subface_values()
506  virtual void
508  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
509  const unsigned int face_no,
510  const unsigned int subface_no,
511  const Quadrature<dim - 1> &quadrature,
512  const typename Mapping<dim, spacedim>::InternalDataBase &internal_data,
514  &output_data) const override;
515 
516  // documentation can be found in Mapping::fill_fe_immersed_surface_values()
517  virtual void
519  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
521  const typename Mapping<dim, spacedim>::InternalDataBase &internal_data,
523  &output_data) const override;
524 
533  const unsigned int polynomial_degree;
534 
535  /*
536  * The default line support points. These are used when computing the
537  * location in real space of the support points on lines and quads, which
538  * are needed by the Manifold<dim,spacedim> class.
539  *
540  * The number of points depends on the degree of this class, and it matches
541  * the number of degrees of freedom of an FE_Q<1>(this->degree).
542  */
543  const std::vector<Point<1>> line_support_points;
544 
545  /*
546  * The one-dimensional polynomials defined as Lagrange polynomials from the
547  * line support points. These are used for point evaluations and match the
548  * polynomial space of an FE_Q<1>(this->degree).
549  */
550  const std::vector<Polynomials::Polynomial<double>> polynomials_1d;
551 
552  /*
553  * The numbering from the lexicographic to the hierarchical ordering used
554  * when expanding the tensor product with the mapping support points (which
555  * come in hierarchical numbers).
556  */
557  const std::vector<unsigned int> renumber_lexicographic_to_hierarchic;
558 
559  /*
560  * The support points in reference coordinates. These are used for
561  * constructing approximations of the output of
562  * compute_mapping_support_points() when evaluating the mapping on the fly,
563  * rather than going through the FEValues interface provided by
564  * InternalData.
565  *
566  * The number of points depends on the degree of this class, and it matches
567  * the number of degrees of freedom of an FE_Q<dim>(this->degree).
568  */
569  const std::vector<Point<dim>> unit_cell_support_points;
570 
590  const std::vector<Table<2, double>>
592 
606 
636  virtual std::vector<Point<spacedim>>
638  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const;
639 
644  Point<dim>
646  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
647  const Point<spacedim> &p,
648  const Point<dim> &initial_p_unit) const;
649 
663  virtual void
665  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
666  std::vector<Point<spacedim>> &a) const;
667 
682  virtual void
684  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
685  std::vector<Point<spacedim>> &a) const;
686 
687  // Make MappingQCache a friend since it needs to call the
688  // compute_mapping_support_points() function.
689  template <int, int>
690  friend class MappingQCache;
691 };
692 
693 
694 
700 template <int dim, int spacedim = dim>
702 
706 /*----------------------------------------------------------------------*/
707 
708 #ifndef DOXYGEN
709 
710 template <int dim, int spacedim>
711 inline bool
713 {
714  return true;
715 }
716 
717 #endif // DOXYGEN
718 
719 /* -------------- declaration of explicit specializations ------------- */
720 
721 
723 
724 #endif
Triangulation< dim, spacedim >::cell_iterator cell_of_current_support_points
Definition: mapping_q.h:445
QGaussLobatto< 1 > line_support_points
Definition: mapping_q.h:391
virtual std::size_t memory_consumption() const override
Definition: mapping_q.cc:64
std::vector< Point< spacedim > > mapping_support_points
Definition: mapping_q.h:439
const unsigned int polynomial_degree
Definition: mapping_q.h:370
std::array< std::vector< Tensor< 1, dim > >, GeometryInfo< dim >::faces_per_cell *(dim - 1)> unit_tangentials
Definition: mapping_q.h:364
internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > * output_data
Definition: mapping_q.h:459
std::vector< AlignedVector< Tensor< 1, spacedim > > > aux
Definition: mapping_q.h:433
void initialize_face(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
Definition: mapping_q.cc:162
InternalData(const unsigned int polynomial_degree)
Definition: mapping_q.cc:51
AlignedVector< double > volume_elements
Definition: mapping_q.h:451
void initialize(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
Definition: mapping_q.cc:82
std::vector< Point< dim > > quadrature_points
Definition: mapping_q.h:347
const unsigned int n_shape_functions
Definition: mapping_q.h:381
internal::MatrixFreeFunctions::ShapeInfo< VectorizedArrayType > shape_info
Definition: mapping_q.h:417
AlignedVector< VectorizedArrayType > scratch
Definition: mapping_q.h:423
virtual void add_line_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim >> &a) const
Definition: mapping_q.cc:1458
const std::vector< unsigned int > renumber_lexicographic_to_hierarchic
Definition: mapping_q.h:557
const Table< 2, double > support_point_weights_cell
Definition: mapping_q.h:605
virtual std::vector< Point< spacedim > > compute_mapping_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const
Definition: mapping_q.cc:1635
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const override
Definition: mapping_q.cc:731
void fill_mapping_data_for_generic_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const ArrayView< const Point< dim >> &unit_points, const UpdateFlags update_flags, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const
Definition: mapping_q.cc:1241
void fill_mapping_data_for_face_quadrature(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_number, const Quadrature< dim - 1 > &face_quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const
Definition: mapping_q.cc:1281
virtual BoundingBox< spacedim > get_bounding_box(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const override
Definition: mapping_q.cc:1734
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags flags, const Quadrature< dim - 1 > &quadrature) const override
Definition: mapping_q.cc:767
virtual void fill_fe_subface_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int subface_no, const Quadrature< dim - 1 > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
Definition: mapping_q.cc:1057
virtual CellSimilarity::Similarity fill_fe_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const Quadrature< dim > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
Definition: mapping_q.cc:786
virtual void fill_fe_immersed_surface_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const NonMatching::ImmersedSurfaceQuadrature< dim > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
Definition: mapping_q.cc:1110
const unsigned int polynomial_degree
Definition: mapping_q.h:533
virtual bool preserves_vertex_locations() const override
virtual void transform(const ArrayView< const Tensor< 1, dim >> &input, const MappingKind kind, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const ArrayView< Tensor< 1, spacedim >> &output) const override
Definition: mapping_q.cc:1317
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_face_data(const UpdateFlags flags, const hp::QCollection< dim - 1 > &quadrature) const override
Definition: mapping_q.cc:746
Point< dim > transform_real_to_unit_cell_internal(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p, const Point< dim > &initial_p_unit) const
Definition: mapping_q.cc:324
const std::vector< Table< 2, double > > support_point_weights_perimeter_to_interior
Definition: mapping_q.h:591
const std::vector< Point< 1 > > line_support_points
Definition: mapping_q.h:543
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const override
Definition: mapping_q.cc:290
const std::vector< Polynomials::Polynomial< double > > polynomials_1d
Definition: mapping_q.h:550
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
Definition: mapping_q.cc:272
const std::vector< Point< dim > > unit_cell_support_points
Definition: mapping_q.h:569
virtual void transform_points_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const ArrayView< const Point< spacedim >> &real_points, const ArrayView< Point< dim >> &unit_points) const override
Definition: mapping_q.cc:593
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const override
Definition: mapping_q.cc:473
MappingQ(const unsigned int polynomial_degree)
Definition: mapping_q.cc:220
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
Definition: mapping_q.cc:675
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const hp::QCollection< dim - 1 > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
Definition: mapping_q.cc:1005
virtual bool is_compatible_with(const ReferenceCell &reference_cell) const override
Definition: mapping_q.cc:1744
virtual void add_quad_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim >> &a) const
Definition: mapping_q.cc:1624
unsigned int get_degree() const
Definition: mapping_q.cc:281
Abstract base class for mapping classes.
Definition: mapping.h:317
#define DEAL_II_DEPRECATED
Definition: config.h:178
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
UpdateFlags
MappingKind
Definition: mapping.h:78
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)