Reference documentation for deal.II version 9.2.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 - 2020 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 
24 #include <deal.II/base/table.h>
26 
27 #include <deal.II/fe/mapping.h>
28 
30 
32 
33 #include <array>
34 #include <cmath>
35 
37 
38 template <int, int>
39 class MappingQ;
40 
41 template <int, int>
42 class MappingQCache;
43 
44 
47 
48 
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 
173  // for documentation, see the Mapping base class
174  virtual Point<spacedim>
176  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
177  const Point<dim> &p) const override;
178 
179  // for documentation, see the Mapping base class
180  virtual Point<dim>
182  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
183  const Point<spacedim> &p) const override;
184 
194  // for documentation, see the Mapping base class
195  virtual void
196  transform(const ArrayView<const Tensor<1, dim>> & input,
197  const MappingKind kind,
198  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
199  const ArrayView<Tensor<1, spacedim>> &output) const override;
200 
201  // for documentation, see the Mapping base class
202  virtual void
204  const MappingKind kind,
205  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
206  const ArrayView<Tensor<2, spacedim>> &output) const override;
207 
208  // for documentation, see the Mapping base class
209  virtual void
210  transform(const ArrayView<const Tensor<2, dim>> & input,
211  const MappingKind kind,
212  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
213  const ArrayView<Tensor<2, spacedim>> &output) const override;
214 
215  // for documentation, see the Mapping base class
216  virtual void
218  const MappingKind kind,
219  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
220  const ArrayView<Tensor<3, spacedim>> &output) const override;
221 
222  // for documentation, see the Mapping base class
223  virtual void
224  transform(const ArrayView<const Tensor<3, dim>> & input,
225  const MappingKind kind,
226  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
227  const ArrayView<Tensor<3, spacedim>> &output) const override;
228 
238 public:
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 
300 
305  const double &
306  shape(const unsigned int qpoint, const unsigned int shape_nr) const;
307 
311  double &
312  shape(const unsigned int qpoint, const unsigned int shape_nr);
313 
317  const Tensor<1, dim> &
318  derivative(const unsigned int qpoint, const unsigned int shape_nr) const;
319 
324  derivative(const unsigned int qpoint, const unsigned int shape_nr);
325 
329  const Tensor<2, dim> &
330  second_derivative(const unsigned int qpoint,
331  const unsigned int shape_nr) const;
332 
337  second_derivative(const unsigned int qpoint, const unsigned int shape_nr);
338 
342  const Tensor<3, dim> &
343  third_derivative(const unsigned int qpoint,
344  const unsigned int shape_nr) const;
345 
350  third_derivative(const unsigned int qpoint, const unsigned int shape_nr);
351 
355  const Tensor<4, dim> &
356  fourth_derivative(const unsigned int qpoint,
357  const unsigned int shape_nr) const;
358 
363  fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr);
364 
368  virtual std::size_t
369  memory_consumption() const override;
370 
376  std::vector<double> shape_values;
377 
383  std::vector<Tensor<1, dim>> shape_derivatives;
384 
391  std::vector<Tensor<2, dim>> shape_second_derivatives;
392 
399  std::vector<Tensor<3, dim>> shape_third_derivatives;
400 
407  std::vector<Tensor<4, dim>> shape_fourth_derivatives;
408 
422  std::array<std::vector<Tensor<1, dim>>,
425 
430  const unsigned int polynomial_degree;
431 
441  const unsigned int n_shape_functions;
442 
443  /*
444  * The default line support points. Is used in when the shape function
445  * values are computed.
446  *
447  * The number of quadrature points depends on the degree of this
448  * class, and it matches the number of degrees of freedom of an
449  * FE_Q<1>(this->degree).
450  */
452 
460 
466 
472 
478 
484 
490 
495 
505  mutable std::vector<DerivativeForm<1, dim, spacedim>> covariant;
506 
514  mutable std::vector<DerivativeForm<1, dim, spacedim>> contravariant;
515 
519  mutable std::vector<std::vector<Tensor<1, spacedim>>> aux;
520 
525  mutable std::vector<Point<spacedim>> mapping_support_points;
526 
532 
537  mutable std::vector<double> volume_elements;
538  };
539 
540 
541  // documentation can be found in Mapping::requires_update_flags()
542  virtual UpdateFlags
543  requires_update_flags(const UpdateFlags update_flags) const override;
544 
545  // documentation can be found in Mapping::get_data()
546  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
547  get_data(const UpdateFlags, const Quadrature<dim> &quadrature) const override;
548 
549  // documentation can be found in Mapping::get_face_data()
550  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
551  get_face_data(const UpdateFlags flags,
552  const Quadrature<dim - 1> &quadrature) const override;
553 
554  // documentation can be found in Mapping::get_subface_data()
555  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
556  get_subface_data(const UpdateFlags flags,
557  const Quadrature<dim - 1> &quadrature) const override;
558 
559  // documentation can be found in Mapping::fill_fe_values()
562  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
563  const CellSimilarity::Similarity cell_similarity,
564  const Quadrature<dim> & quadrature,
565  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
567  &output_data) const override;
568 
569  // documentation can be found in Mapping::fill_fe_face_values()
570  virtual void
572  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
573  const unsigned int face_no,
574  const Quadrature<dim - 1> & quadrature,
575  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
577  &output_data) const override;
578 
579  // documentation can be found in Mapping::fill_fe_subface_values()
580  virtual void
582  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
583  const unsigned int face_no,
584  const unsigned int subface_no,
585  const Quadrature<dim - 1> & quadrature,
586  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
588  &output_data) const override;
589 
594 protected:
599  const unsigned int polynomial_degree;
600 
601  /*
602  * The default line support points. These are used when computing
603  * the location in real space of the support points on lines and
604  * quads, which are asked to the Manifold<dim,spacedim> class.
605  *
606  * The number of quadrature points depends on the degree of this
607  * class, and it matches the number of degrees of freedom of an
608  * FE_Q<1>(this->degree).
609  */
611 
631  std::vector<Table<2, double>> support_point_weights_perimeter_to_interior;
632 
646 
676  virtual std::vector<Point<spacedim>>
678  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const;
679 
684  Point<dim>
686  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
687  const Point<spacedim> & p,
688  const Point<dim> &initial_p_unit) const;
689 
703  virtual void
705  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
706  std::vector<Point<spacedim>> & a) const;
707 
722  virtual void
724  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
725  std::vector<Point<spacedim>> & a) const;
726 
727  // Make MappingQ a friend since it needs to call the fill_fe_values()
728  // functions on its MappingQGeneric(1) sub-object.
729  template <int, int>
730  friend class MappingQ;
731 
732  // Make MappingQCache a friend since it needs to call the
733  // compute_mapping_support_points() function.
734  template <int, int>
735  friend class MappingQCache;
736 };
737 
738 
739 
742 /*----------------------------------------------------------------------*/
743 
744 #ifndef DOXYGEN
745 
746 template <int dim, int spacedim>
747 inline const double &
749  const unsigned int qpoint,
750  const unsigned int shape_nr) const
751 {
752  AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
753  return shape_values[qpoint * n_shape_functions + shape_nr];
754 }
755 
756 
757 
758 template <int dim, int spacedim>
759 inline double &
761  const unsigned int shape_nr)
762 {
763  AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
764  return shape_values[qpoint * n_shape_functions + shape_nr];
765 }
766 
767 
768 template <int dim, int spacedim>
769 inline const Tensor<1, dim> &
771  const unsigned int qpoint,
772  const unsigned int shape_nr) const
773 {
774  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
775  shape_derivatives.size());
776  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
777 }
778 
779 
780 
781 template <int dim, int spacedim>
782 inline Tensor<1, dim> &
784  const unsigned int qpoint,
785  const unsigned int shape_nr)
786 {
787  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
788  shape_derivatives.size());
789  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
790 }
791 
792 
793 template <int dim, int spacedim>
794 inline const Tensor<2, dim> &
796  const unsigned int qpoint,
797  const unsigned int shape_nr) const
798 {
799  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
800  shape_second_derivatives.size());
801  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
802 }
803 
804 
805 template <int dim, int spacedim>
806 inline Tensor<2, dim> &
808  const unsigned int qpoint,
809  const unsigned int shape_nr)
810 {
811  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
812  shape_second_derivatives.size());
813  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
814 }
815 
816 template <int dim, int spacedim>
817 inline const Tensor<3, dim> &
819  const unsigned int qpoint,
820  const unsigned int shape_nr) const
821 {
822  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
823  shape_third_derivatives.size());
824  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
825 }
826 
827 
828 template <int dim, int spacedim>
829 inline Tensor<3, dim> &
831  const unsigned int qpoint,
832  const unsigned int shape_nr)
833 {
834  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
835  shape_third_derivatives.size());
836  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
837 }
838 
839 
840 template <int dim, int spacedim>
841 inline const Tensor<4, dim> &
843  const unsigned int qpoint,
844  const unsigned int shape_nr) const
845 {
846  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
847  shape_fourth_derivatives.size());
848  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
849 }
850 
851 
852 template <int dim, int spacedim>
853 inline Tensor<4, dim> &
855  const unsigned int qpoint,
856  const unsigned int shape_nr)
857 {
858  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
859  shape_fourth_derivatives.size());
860  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
861 }
862 
863 
864 
865 template <int dim, int spacedim>
866 inline bool
868 {
869  return true;
870 }
871 
872 #endif // DOXYGEN
873 
874 /* -------------- declaration of explicit specializations ------------- */
875 
876 
878 
879 #endif
std::vector< Tensor< 2, dim > > shape_second_derivatives
const unsigned int polynomial_degree
Table< 2, double > support_point_weights_cell
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:1649
QGaussLobatto< 1 > line_support_points
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
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)
AlignedVector< VectorizedArray< double > > scratch
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< dim - 1 > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, ::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
MappingKind
Definition: mapping.h:62
InternalData(const unsigned int polynomial_degree)
const unsigned int polynomial_degree
QGaussLobatto< 1 > line_support_points
std::vector< Table< 2, double > > support_point_weights_perimeter_to_interior
void compute_shape_function_values(const std::vector< Point< dim >> &unit_points)
std::vector< std::vector< Tensor< 1, spacedim > > > aux
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:302
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:359
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
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
AlignedVector< VectorizedArray< double > > hessians_quad
AlignedVector< VectorizedArray< double > > values_dofs
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
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:358
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
std::vector< DerivativeForm< 1, dim, spacedim > > contravariant
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
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_face_data(const UpdateFlags flags, const Quadrature< dim - 1 > &quadrature) const override