Reference documentation for deal.II version Git 8ba89c874e 2021-08-04 22:06:33 -0600
\(\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 - 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 #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 
34 
35 #include <array>
36 #include <cmath>
37 
39 
40 template <int, int>
41 class MappingQ;
42 
43 template <int, int>
45 
46 
49 
50 
111 template <int dim, int spacedim = dim>
112 class MappingQ : public Mapping<dim, spacedim>
113 {
114 public:
120  MappingQ(const unsigned int polynomial_degree);
121 
127  DEAL_II_DEPRECATED_EARLY
128  MappingQ(const unsigned int polynomial_degree,
129  const bool use_mapping_q_on_all_cells);
130 
134  MappingQ(const MappingQ<dim, spacedim> &mapping);
135 
136  // for documentation, see the Mapping base class
137  virtual std::unique_ptr<Mapping<dim, spacedim>>
138  clone() const override;
139 
144  unsigned int
145  get_degree() const;
146 
151  virtual bool
152  preserves_vertex_locations() const override;
153 
154  // for documentation, see the Mapping base class
155  virtual BoundingBox<spacedim>
157  &cell) const override;
158 
159  virtual bool
160  is_compatible_with(const ReferenceCell &reference_cell) const override;
161 
167  // for documentation, see the Mapping base class
168  virtual Point<spacedim>
170  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
171  const Point<dim> &p) const override;
172 
173  // for documentation, see the Mapping base class
174  virtual Point<dim>
176  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
177  const Point<spacedim> &p) const override;
178 
179  // for documentation, see the Mapping base class
180  virtual void
182  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
183  const ArrayView<const Point<spacedim>> & real_points,
184  const ArrayView<Point<dim>> &unit_points) const override;
185 
195  // for documentation, see the Mapping base class
196  virtual void
197  transform(const ArrayView<const Tensor<1, dim>> & input,
198  const MappingKind kind,
199  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
200  const ArrayView<Tensor<1, spacedim>> &output) const override;
201 
202  // for documentation, see the Mapping base class
203  virtual void
205  const MappingKind kind,
206  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
207  const ArrayView<Tensor<2, spacedim>> &output) const override;
208 
209  // for documentation, see the Mapping base class
210  virtual void
211  transform(const ArrayView<const Tensor<2, dim>> & input,
212  const MappingKind kind,
213  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
214  const ArrayView<Tensor<2, spacedim>> &output) const override;
215 
216  // for documentation, see the Mapping base class
217  virtual void
219  const MappingKind kind,
220  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
221  const ArrayView<Tensor<3, spacedim>> &output) const override;
222 
223  // for documentation, see the Mapping base class
224  virtual void
225  transform(const ArrayView<const Tensor<3, dim>> & input,
226  const MappingKind kind,
227  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
228  const ArrayView<Tensor<3, spacedim>> &output) const override;
229 
250  class InternalData : public Mapping<dim, spacedim>::InternalDataBase
251  {
252  public:
257  InternalData(const unsigned int polynomial_degree);
258 
267  void
268  initialize(const UpdateFlags update_flags,
269  const Quadrature<dim> &quadrature,
270  const unsigned int n_original_q_points);
271 
277  void
278  initialize_face(const UpdateFlags update_flags,
279  const Quadrature<dim> &quadrature,
280  const unsigned int n_original_q_points);
281 
297  void
298  compute_shape_function_values(const std::vector<Point<dim>> &unit_points);
299 
304  const double &
305  shape(const unsigned int qpoint, const unsigned int shape_nr) const;
306 
310  double &
311  shape(const unsigned int qpoint, const unsigned int shape_nr);
312 
316  const Tensor<1, dim> &
317  derivative(const unsigned int qpoint, const unsigned int shape_nr) const;
318 
323  derivative(const unsigned int qpoint, const unsigned int shape_nr);
324 
328  const Tensor<2, dim> &
329  second_derivative(const unsigned int qpoint,
330  const unsigned int shape_nr) const;
331 
336  second_derivative(const unsigned int qpoint, const unsigned int shape_nr);
337 
341  const Tensor<3, dim> &
342  third_derivative(const unsigned int qpoint,
343  const unsigned int shape_nr) const;
344 
349  third_derivative(const unsigned int qpoint, const unsigned int shape_nr);
350 
354  const Tensor<4, dim> &
355  fourth_derivative(const unsigned int qpoint,
356  const unsigned int shape_nr) const;
357 
362  fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr);
363 
367  virtual std::size_t
368  memory_consumption() const override;
369 
376 
383 
391 
399 
407 
421  std::array<std::vector<Tensor<1, dim>>,
424 
429  const unsigned int polynomial_degree;
430 
440  const unsigned int n_shape_functions;
441 
442  /*
443  * The default line support points. Is used in when the shape function
444  * values are computed.
445  *
446  * The number of quadrature points depends on the degree of this
447  * class, and it matches the number of degrees of freedom of an
448  * FE_Q<1>(this->degree).
449  */
451 
459 
465 
471 
477 
483 
489 
494 
505 
514 
518  mutable std::vector<AlignedVector<Tensor<1, spacedim>>> aux;
519 
524  mutable std::vector<Point<spacedim>> mapping_support_points;
525 
531 
537  };
538 
539 
540  // documentation can be found in Mapping::requires_update_flags()
541  virtual UpdateFlags
542  requires_update_flags(const UpdateFlags update_flags) const override;
543 
544  // documentation can be found in Mapping::get_data()
545  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
546  get_data(const UpdateFlags, const Quadrature<dim> &quadrature) const override;
547 
549 
550  // documentation can be found in Mapping::get_face_data()
551  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
552  get_face_data(const UpdateFlags flags,
553  const hp::QCollection<dim - 1> &quadrature) const override;
554 
555  // documentation can be found in Mapping::get_subface_data()
556  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
557  get_subface_data(const UpdateFlags flags,
558  const Quadrature<dim - 1> &quadrature) const override;
559 
560  // documentation can be found in Mapping::fill_fe_values()
563  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
564  const CellSimilarity::Similarity cell_similarity,
565  const Quadrature<dim> & quadrature,
566  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
568  &output_data) const override;
569 
571 
572  // documentation can be found in Mapping::fill_fe_face_values()
573  virtual void
575  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
576  const unsigned int face_no,
577  const hp::QCollection<dim - 1> & quadrature,
578  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
580  &output_data) const override;
581 
582  // documentation can be found in Mapping::fill_fe_subface_values()
583  virtual void
585  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
586  const unsigned int face_no,
587  const unsigned int subface_no,
588  const Quadrature<dim - 1> & quadrature,
589  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
591  &output_data) const override;
592 
593 
611  void
613  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
614  const ArrayView<const Point<dim>> & unit_points,
615  const UpdateFlags update_flags,
617  &output_data) const;
618 
623 protected:
628  const unsigned int polynomial_degree;
629 
630  /*
631  * The default line support points. These are used when computing the
632  * location in real space of the support points on lines and quads, which
633  * are needed by the Manifold<dim,spacedim> class.
634  *
635  * The number of points depends on the degree of this class, and it matches
636  * the number of degrees of freedom of an FE_Q<1>(this->degree).
637  */
638  const std::vector<Point<1>> line_support_points;
639 
640  /*
641  * The one-dimensional polynomials defined as Lagrange polynomials from the
642  * line support points. These are used for point evaluations and match the
643  * polynomial space of an FE_Q<1>(this->degree).
644  */
645  const std::vector<Polynomials::Polynomial<double>> polynomials_1d;
646 
647  /*
648  * The numbering from the lexicographic to the hierarchical ordering used
649  * when expanding the tensor product with the mapping support points (which
650  * come in hierarchical numbers).
651  */
652  const std::vector<unsigned int> renumber_lexicographic_to_hierarchic;
653 
654  /*
655  * The support points in reference coordinates. These are used for
656  * constructing approximations of the output of
657  * compute_mapping_support_points() when evaluating the mapping on the fly,
658  * rather than going through the FEValues interface provided by
659  * InternalData.
660  *
661  * The number of points depends on the degree of this class, and it matches
662  * the number of degrees of freedom of an FE_Q<dim>(this->degree).
663  */
664  const std::vector<Point<dim>> unit_cell_support_points;
665 
685  const std::vector<Table<2, double>>
687 
701 
731  virtual std::vector<Point<spacedim>>
733  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const;
734 
739  Point<dim>
741  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
742  const Point<spacedim> & p,
743  const Point<dim> &initial_p_unit) const;
744 
758  virtual void
760  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
761  std::vector<Point<spacedim>> & a) const;
762 
777  virtual void
779  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
780  std::vector<Point<spacedim>> & a) const;
781 
782  // Make MappingQCache a friend since it needs to call the
783  // compute_mapping_support_points() function.
784  template <int, int>
785  friend class MappingQCache;
786 };
787 
788 
789 
795 template <int dim, int spacedim = dim>
797 
801 /*----------------------------------------------------------------------*/
802 
803 #ifndef DOXYGEN
804 
805 template <int dim, int spacedim>
806 inline const double &
807 MappingQ<dim, spacedim>::InternalData::shape(const unsigned int qpoint,
808  const unsigned int shape_nr) const
809 {
810  AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
811  return shape_values[qpoint * n_shape_functions + shape_nr];
812 }
813 
814 
815 
816 template <int dim, int spacedim>
817 inline double &
818 MappingQ<dim, spacedim>::InternalData::shape(const unsigned int qpoint,
819  const unsigned int shape_nr)
820 {
821  AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
822  return shape_values[qpoint * n_shape_functions + shape_nr];
823 }
824 
825 
826 template <int dim, int spacedim>
827 inline const Tensor<1, dim> &
829  const unsigned int qpoint,
830  const unsigned int shape_nr) const
831 {
832  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
834  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
835 }
836 
837 
838 
839 template <int dim, int spacedim>
840 inline Tensor<1, dim> &
841 MappingQ<dim, spacedim>::InternalData::derivative(const unsigned int qpoint,
842  const unsigned int shape_nr)
843 {
844  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
846  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
847 }
848 
849 
850 template <int dim, int spacedim>
851 inline const Tensor<2, dim> &
853  const unsigned int qpoint,
854  const unsigned int shape_nr) const
855 {
856  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
858  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
859 }
860 
861 
862 template <int dim, int spacedim>
863 inline Tensor<2, dim> &
865  const unsigned int qpoint,
866  const unsigned int shape_nr)
867 {
868  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
870  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
871 }
872 
873 template <int dim, int spacedim>
874 inline const Tensor<3, dim> &
876  const unsigned int qpoint,
877  const unsigned int shape_nr) const
878 {
879  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
881  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
882 }
883 
884 
885 template <int dim, int spacedim>
886 inline Tensor<3, dim> &
888  const unsigned int qpoint,
889  const unsigned int shape_nr)
890 {
891  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
893  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
894 }
895 
896 
897 template <int dim, int spacedim>
898 inline const Tensor<4, dim> &
900  const unsigned int qpoint,
901  const unsigned int shape_nr) const
902 {
903  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
905  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
906 }
907 
908 
909 template <int dim, int spacedim>
910 inline Tensor<4, dim> &
912  const unsigned int qpoint,
913  const unsigned int shape_nr)
914 {
915  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
917  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
918 }
919 
920 
921 
922 template <int dim, int spacedim>
923 inline bool
925 {
926  return true;
927 }
928 
929 #endif // DOXYGEN
930 
931 /* -------------- declaration of explicit specializations ------------- */
932 
933 
935 
936 #endif
const unsigned int polynomial_degree
Definition: mapping_q.h:429
unsigned int get_degree() const
Definition: mapping_q.cc:444
const std::vector< Point< dim > > unit_cell_support_points
Definition: mapping_q.h:664
AlignedVector< Tensor< 4, dim > > shape_fourth_derivatives
Definition: mapping_q.h:406
virtual std::size_t memory_consumption() const override
Definition: mapping_q.cc:64
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:752
AlignedVector< VectorizedArray< double > > values_quad
Definition: mapping_q.h:476
void compute_shape_function_values(const std::vector< Point< dim >> &unit_points)
Definition: mapping_q.cc:271
AlignedVector< DerivativeForm< 1, dim, spacedim > > contravariant
Definition: mapping_q.h:513
const double & shape(const unsigned int qpoint, const unsigned int shape_nr) const
MappingQ(const unsigned int polynomial_degree)
Definition: mapping_q.cc:362
InternalData(const unsigned int polynomial_degree)
Definition: mapping_q.cc:52
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const override
Definition: mapping_q.cc:890
std::array< std::vector< Tensor< 1, dim > >, GeometryInfo< dim >::faces_per_cell *(dim - 1)> unit_tangentials
Definition: mapping_q.h:423
#define AssertIndexRange(index, range)
Definition: exceptions.h:1724
AlignedVector< double > volume_elements
Definition: mapping_q.h:536
virtual void add_quad_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim >> &a) const
Definition: mapping_q.cc:1669
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:1365
AlignedVector< VectorizedArray< double > > gradients_quad
Definition: mapping_q.h:482
MappingKind
Definition: mapping.h:64
const std::vector< unsigned int > renumber_lexicographic_to_hierarchic
Definition: mapping_q.h:652
AlignedVector< VectorizedArray< double > > hessians_quad
Definition: mapping_q.h:488
void initialize_face(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
Definition: mapping_q.cc:215
std::vector< AlignedVector< Tensor< 1, spacedim > > > aux
Definition: mapping_q.h:518
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:905
void initialize(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
Definition: mapping_q.cc:85
const Tensor< 3, dim > & third_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
const unsigned int n_shape_functions
Definition: mapping_q.h:440
AlignedVector< double > shape_values
Definition: mapping_q.h:375
const Tensor< 2, dim > & second_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
AlignedVector< VectorizedArray< double > > values_dofs
Definition: mapping_q.h:470
const Tensor< 1, dim > & derivative(const unsigned int qpoint, const unsigned int shape_nr) const
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:1165
UpdateFlags
const Table< 2, double > support_point_weights_cell
Definition: mapping_q.h:700
const std::vector< Point< 1 > > line_support_points
Definition: mapping_q.h:638
Abstract base class for mapping classes.
Definition: mapping.h:303
internal::MatrixFreeFunctions::ShapeInfo< VectorizedArray< double > > shape_info
Definition: mapping_q.h:458
AlignedVector< VectorizedArray< double > > scratch
Definition: mapping_q.h:464
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
const unsigned int polynomial_degree
Definition: mapping_q.h:628
const Tensor< 4, dim > & fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
AlignedVector< Tensor< 2, dim > > shape_second_derivatives
Definition: mapping_q.h:390
virtual BoundingBox< spacedim > get_bounding_box(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const override
Definition: mapping_q.cc:1779
const std::vector< Polynomials::Polynomial< double > > polynomials_1d
Definition: mapping_q.h:645
QGaussLobatto< 1 > line_support_points
Definition: mapping_q.h:450
Triangulation< dim, spacedim >::cell_iterator cell_of_current_support_points
Definition: mapping_q.h:530
std::vector< Point< spacedim > > mapping_support_points
Definition: mapping_q.h:524
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:1214
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
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:926
virtual bool is_compatible_with(const ReferenceCell &reference_cell) const override
Definition: mapping_q.cc:1789
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:1264
virtual std::vector< Point< spacedim > > compute_mapping_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const
Definition: mapping_q.cc:1680
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:453
AlignedVector< Tensor< 3, dim > > shape_third_derivatives
Definition: mapping_q.h:398
virtual bool preserves_vertex_locations() const override
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
AlignedVector< Tensor< 1, dim > > shape_derivatives
Definition: mapping_q.h:382
const std::vector< Table< 2, double > > support_point_weights_perimeter_to_interior
Definition: mapping_q.h:686
size_type size() const
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:488
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
Definition: mapping_q.cc:435
virtual void add_line_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim >> &a) const
Definition: mapping_q.cc:1503
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
Definition: mapping_q.cc:834
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:945
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:631
AlignedVector< DerivativeForm< 1, dim, spacedim > > covariant
Definition: mapping_q.h:504