Reference documentation for deal.II version 9.3.0
\(\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_generic.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_generic_h
17 #define dealii_mapping_q_generic_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>
44 class MappingQCache;
45 
46 
49 
50 
134 template <int dim, int spacedim = dim>
135 class MappingQGeneric : public Mapping<dim, spacedim>
136 {
137 public:
143  MappingQGeneric(const unsigned int polynomial_degree);
144 
149 
150  // for documentation, see the Mapping base class
151  virtual std::unique_ptr<Mapping<dim, spacedim>>
152  clone() const override;
153 
158  unsigned int
159  get_degree() const;
160 
165  virtual bool
166  preserves_vertex_locations() const override;
167 
168  // for documentation, see the Mapping base class
169  virtual BoundingBox<spacedim>
171  &cell) const override;
172 
173  virtual bool
174  is_compatible_with(const ReferenceCell &reference_cell) const override;
175 
181  // for documentation, see the Mapping base class
182  virtual Point<spacedim>
184  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
185  const Point<dim> &p) const override;
186 
187  // for documentation, see the Mapping base class
188  virtual Point<dim>
190  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
191  const Point<spacedim> &p) const override;
192 
193  // for documentation, see the Mapping base class
194  virtual void
196  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
197  const ArrayView<const Point<spacedim>> & real_points,
198  const ArrayView<Point<dim>> &unit_points) const override;
199 
209  // for documentation, see the Mapping base class
210  virtual void
211  transform(const ArrayView<const Tensor<1, dim>> & input,
212  const MappingKind kind,
213  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
214  const ArrayView<Tensor<1, 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<2, spacedim>> &output) const override;
222 
223  // for documentation, see the Mapping base class
224  virtual void
225  transform(const ArrayView<const Tensor<2, dim>> & input,
226  const MappingKind kind,
227  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
228  const ArrayView<Tensor<2, spacedim>> &output) const override;
229 
230  // for documentation, see the Mapping base class
231  virtual void
233  const MappingKind kind,
234  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
235  const ArrayView<Tensor<3, spacedim>> &output) const override;
236 
237  // for documentation, see the Mapping base class
238  virtual void
239  transform(const ArrayView<const Tensor<3, dim>> & input,
240  const MappingKind kind,
241  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
242  const ArrayView<Tensor<3, spacedim>> &output) const override;
243 
264  class InternalData : public Mapping<dim, spacedim>::InternalDataBase
265  {
266  public:
271  InternalData(const unsigned int polynomial_degree);
272 
281  void
282  initialize(const UpdateFlags update_flags,
283  const Quadrature<dim> &quadrature,
284  const unsigned int n_original_q_points);
285 
291  void
292  initialize_face(const UpdateFlags update_flags,
293  const Quadrature<dim> &quadrature,
294  const unsigned int n_original_q_points);
295 
311  void
312  compute_shape_function_values(const std::vector<Point<dim>> &unit_points);
313 
318  const double &
319  shape(const unsigned int qpoint, const unsigned int shape_nr) const;
320 
324  double &
325  shape(const unsigned int qpoint, const unsigned int shape_nr);
326 
330  const Tensor<1, dim> &
331  derivative(const unsigned int qpoint, const unsigned int shape_nr) const;
332 
337  derivative(const unsigned int qpoint, const unsigned int shape_nr);
338 
342  const Tensor<2, dim> &
343  second_derivative(const unsigned int qpoint,
344  const unsigned int shape_nr) const;
345 
350  second_derivative(const unsigned int qpoint, const unsigned int shape_nr);
351 
355  const Tensor<3, dim> &
356  third_derivative(const unsigned int qpoint,
357  const unsigned int shape_nr) const;
358 
363  third_derivative(const unsigned int qpoint, const unsigned int shape_nr);
364 
368  const Tensor<4, dim> &
369  fourth_derivative(const unsigned int qpoint,
370  const unsigned int shape_nr) const;
371 
376  fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr);
377 
381  virtual std::size_t
382  memory_consumption() const override;
383 
389  std::vector<double> shape_values;
390 
396  std::vector<Tensor<1, dim>> shape_derivatives;
397 
404  std::vector<Tensor<2, dim>> shape_second_derivatives;
405 
412  std::vector<Tensor<3, dim>> shape_third_derivatives;
413 
420  std::vector<Tensor<4, dim>> shape_fourth_derivatives;
421 
435  std::array<std::vector<Tensor<1, dim>>,
438 
443  const unsigned int polynomial_degree;
444 
454  const unsigned int n_shape_functions;
455 
456  /*
457  * The default line support points. Is used in when the shape function
458  * values are computed.
459  *
460  * The number of quadrature points depends on the degree of this
461  * class, and it matches the number of degrees of freedom of an
462  * FE_Q<1>(this->degree).
463  */
465 
473 
479 
485 
491 
497 
503 
508 
518  mutable std::vector<DerivativeForm<1, dim, spacedim>> covariant;
519 
527  mutable std::vector<DerivativeForm<1, dim, spacedim>> contravariant;
528 
532  mutable std::vector<std::vector<Tensor<1, spacedim>>> aux;
533 
538  mutable std::vector<Point<spacedim>> mapping_support_points;
539 
545 
550  mutable std::vector<double> volume_elements;
551  };
552 
553 
554  // documentation can be found in Mapping::requires_update_flags()
555  virtual UpdateFlags
556  requires_update_flags(const UpdateFlags update_flags) const override;
557 
558  // documentation can be found in Mapping::get_data()
559  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
560  get_data(const UpdateFlags, const Quadrature<dim> &quadrature) const override;
561 
563 
564  // documentation can be found in Mapping::get_face_data()
565  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
566  get_face_data(const UpdateFlags flags,
567  const hp::QCollection<dim - 1> &quadrature) const override;
568 
569  // documentation can be found in Mapping::get_subface_data()
570  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
571  get_subface_data(const UpdateFlags flags,
572  const Quadrature<dim - 1> &quadrature) const override;
573 
574  // documentation can be found in Mapping::fill_fe_values()
577  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
578  const CellSimilarity::Similarity cell_similarity,
579  const Quadrature<dim> & quadrature,
580  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
582  &output_data) const override;
583 
585 
586  // documentation can be found in Mapping::fill_fe_face_values()
587  virtual void
589  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
590  const unsigned int face_no,
591  const hp::QCollection<dim - 1> & quadrature,
592  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
594  &output_data) const override;
595 
596  // documentation can be found in Mapping::fill_fe_subface_values()
597  virtual void
599  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
600  const unsigned int face_no,
601  const unsigned int subface_no,
602  const Quadrature<dim - 1> & quadrature,
603  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
605  &output_data) const override;
606 
607 
625  void
627  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
628  const ArrayView<const Point<dim>> & unit_points,
629  const UpdateFlags update_flags,
631  &output_data) const;
632 
637 protected:
642  const unsigned int polynomial_degree;
643 
644  /*
645  * The default line support points. These are used when computing the
646  * location in real space of the support points on lines and quads, which
647  * are needed by the Manifold<dim,spacedim> class.
648  *
649  * The number of points depends on the degree of this class, and it matches
650  * the number of degrees of freedom of an FE_Q<1>(this->degree).
651  */
652  const std::vector<Point<1>> line_support_points;
653 
654  /*
655  * The one-dimensional polynomials defined as Lagrange polynomials from the
656  * line support points. These are used for point evaluations and match the
657  * polynomial space of an FE_Q<1>(this->degree).
658  */
659  const std::vector<Polynomials::Polynomial<double>> polynomials_1d;
660 
661  /*
662  * The numbering from the lexicographic to the hierarchical ordering used
663  * when expanding the tensor product with the mapping support points (which
664  * come in hierarchical numbers).
665  */
666  const std::vector<unsigned int> renumber_lexicographic_to_hierarchic;
667 
668  /*
669  * The support points in reference coordinates. These are used for
670  * constructing approximations of the output of
671  * compute_mapping_support_points() when evaluating the mapping on the fly,
672  * rather than going through the FEValues interface provided by
673  * InternalData.
674  *
675  * The number of points depends on the degree of this class, and it matches
676  * the number of degrees of freedom of an FE_Q<dim>(this->degree).
677  */
678  const std::vector<Point<dim>> unit_cell_support_points;
679 
699  const std::vector<Table<2, double>>
701 
715 
745  virtual std::vector<Point<spacedim>>
747  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const;
748 
753  Point<dim>
755  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
756  const Point<spacedim> & p,
757  const Point<dim> &initial_p_unit) const;
758 
772  virtual void
774  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
775  std::vector<Point<spacedim>> & a) const;
776 
791  virtual void
793  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
794  std::vector<Point<spacedim>> & a) const;
795 
796  // Make MappingQ a friend since it needs to call the fill_fe_values()
797  // functions on its MappingQGeneric(1) sub-object.
798  template <int, int>
799  friend class MappingQ;
800 
801  // Make MappingQCache a friend since it needs to call the
802  // compute_mapping_support_points() function.
803  template <int, int>
804  friend class MappingQCache;
805 };
806 
807 
808 
811 /*----------------------------------------------------------------------*/
812 
813 #ifndef DOXYGEN
814 
815 template <int dim, int spacedim>
816 inline const double &
818  const unsigned int qpoint,
819  const unsigned int shape_nr) const
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 
827 template <int dim, int spacedim>
828 inline double &
830  const unsigned int shape_nr)
831 {
832  AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
833  return shape_values[qpoint * n_shape_functions + shape_nr];
834 }
835 
836 
837 template <int dim, int spacedim>
838 inline const Tensor<1, dim> &
840  const unsigned int qpoint,
841  const unsigned int shape_nr) const
842 {
843  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
844  shape_derivatives.size());
845  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
846 }
847 
848 
849 
850 template <int dim, int spacedim>
851 inline Tensor<1, dim> &
853  const unsigned int qpoint,
854  const unsigned int shape_nr)
855 {
856  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
857  shape_derivatives.size());
858  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
859 }
860 
861 
862 template <int dim, int spacedim>
863 inline const Tensor<2, dim> &
865  const unsigned int qpoint,
866  const unsigned int shape_nr) const
867 {
868  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
869  shape_second_derivatives.size());
870  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
871 }
872 
873 
874 template <int dim, int spacedim>
875 inline Tensor<2, dim> &
877  const unsigned int qpoint,
878  const unsigned int shape_nr)
879 {
880  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
881  shape_second_derivatives.size());
882  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
883 }
884 
885 template <int dim, int spacedim>
886 inline const Tensor<3, dim> &
888  const unsigned int qpoint,
889  const unsigned int shape_nr) const
890 {
891  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
892  shape_third_derivatives.size());
893  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
894 }
895 
896 
897 template <int dim, int spacedim>
898 inline Tensor<3, dim> &
900  const unsigned int qpoint,
901  const unsigned int shape_nr)
902 {
903  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
904  shape_third_derivatives.size());
905  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
906 }
907 
908 
909 template <int dim, int spacedim>
910 inline const Tensor<4, dim> &
912  const unsigned int qpoint,
913  const unsigned int shape_nr) const
914 {
915  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
916  shape_fourth_derivatives.size());
917  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
918 }
919 
920 
921 template <int dim, int spacedim>
922 inline Tensor<4, dim> &
924  const unsigned int qpoint,
925  const unsigned int shape_nr)
926 {
927  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
928  shape_fourth_derivatives.size());
929  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
930 }
931 
932 
933 
934 template <int dim, int spacedim>
935 inline bool
937 {
938  return true;
939 }
940 
941 #endif // DOXYGEN
942 
943 /* -------------- declaration of explicit specializations ------------- */
944 
945 
947 
948 #endif
std::vector< Tensor< 2, dim > > shape_second_derivatives
const unsigned int polynomial_degree
virtual void add_quad_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim >> &a) const
AlignedVector< VectorizedArray< double > > values_quad
const Tensor< 1, dim > & derivative(const unsigned int qpoint, const unsigned int shape_nr) const
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
virtual void add_line_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim >> &a) const
const Tensor< 3, dim > & third_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const override
#define AssertIndexRange(index, range)
Definition: exceptions.h:1690
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
const Table< 2, double > support_point_weights_cell
AlignedVector< VectorizedArray< double > > gradients_quad
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
MappingQGeneric(const unsigned int polynomial_degree)
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
AlignedVector< VectorizedArray< double > > scratch
MappingKind
Definition: mapping.h:64
InternalData(const unsigned int polynomial_degree)
const unsigned int polynomial_degree
const std::vector< Polynomials::Polynomial< double > > polynomials_1d
QGaussLobatto< 1 > line_support_points
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_face_data(const UpdateFlags flags, const hp::QCollection< dim - 1 > &quadrature) const override
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
virtual bool is_compatible_with(const ReferenceCell &reference_cell) const override
void compute_shape_function_values(const std::vector< Point< dim >> &unit_points)
std::vector< std::vector< Tensor< 1, spacedim > > > aux
const std::vector< Point< dim > > unit_cell_support_points
UpdateFlags
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
Abstract base class for mapping classes.
Definition: mapping.h:303
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:395
std::vector< Point< spacedim > > mapping_support_points
std::vector< Tensor< 3, dim > > shape_third_derivatives
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const override
void initialize(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
std::vector< double > volume_elements
const unsigned int n_shape_functions
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags flags, const Quadrature< dim - 1 > &quadrature) const override
std::vector< Tensor< 1, dim > > shape_derivatives
const std::vector< Table< 2, double > > support_point_weights_perimeter_to_interior
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
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
AlignedVector< VectorizedArray< double > > hessians_quad
AlignedVector< VectorizedArray< double > > values_dofs
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
unsigned int get_degree() const
const double & shape(const unsigned int qpoint, const unsigned int shape_nr) const
std::array< std::vector< Tensor< 1, dim > >, GeometryInfo< dim >::faces_per_cell *(dim - 1)> unit_tangentials
const std::vector< Point< 1 > > line_support_points
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:394
virtual bool preserves_vertex_locations() const override
const Tensor< 2, dim > & second_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
Triangulation< dim, spacedim >::cell_iterator cell_of_current_support_points
virtual std::size_t memory_consumption() const override
virtual std::vector< Point< spacedim > > compute_mapping_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const
const std::vector< unsigned int > renumber_lexicographic_to_hierarchic
std::vector< DerivativeForm< 1, dim, spacedim > > contravariant
virtual BoundingBox< spacedim > get_bounding_box(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const override
std::vector< Tensor< 4, dim > > shape_fourth_derivatives
void initialize_face(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
std::vector< double > shape_values
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const override
const Tensor< 4, dim > & fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
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
std::vector< DerivativeForm< 1, dim, spacedim > > covariant
internal::MatrixFreeFunctions::ShapeInfo< VectorizedArray< double > > shape_info