Reference documentation for deal.II version Git 932f7faded 2020-11-28 20:02:43 +0100
\(\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 
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 template <int, int>
40 class MappingQ;
41 
42 template <int, int>
43 class MappingQCache;
44 
45 
48 
49 
133 template <int dim, int spacedim = dim>
134 class MappingQGeneric : public Mapping<dim, spacedim>
135 {
136 public:
142  MappingQGeneric(const unsigned int polynomial_degree);
143 
148 
149  // for documentation, see the Mapping base class
150  virtual std::unique_ptr<Mapping<dim, spacedim>>
151  clone() const override;
152 
157  unsigned int
158  get_degree() const;
159 
164  virtual bool
165  preserves_vertex_locations() const override;
166 
172  // for documentation, see the Mapping base class
173  virtual Point<spacedim>
175  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
176  const Point<dim> &p) const override;
177 
178  // for documentation, see the Mapping base class
179  virtual Point<dim>
181  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
182  const Point<spacedim> &p) const override;
183 
184  // for documentation, see the Mapping base class
185  virtual void
187  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
188  const ArrayView<const Point<spacedim>> & real_points,
189  const ArrayView<Point<dim>> &unit_points) const override;
190 
200  // for documentation, see the Mapping base class
201  virtual void
202  transform(const ArrayView<const Tensor<1, dim>> & input,
203  const MappingKind kind,
204  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
205  const ArrayView<Tensor<1, spacedim>> &output) const override;
206 
207  // for documentation, see the Mapping base class
208  virtual void
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
216  transform(const ArrayView<const Tensor<2, dim>> & input,
217  const MappingKind kind,
218  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
219  const ArrayView<Tensor<2, spacedim>> &output) const override;
220 
221  // for documentation, see the Mapping base class
222  virtual void
224  const MappingKind kind,
225  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
226  const ArrayView<Tensor<3, spacedim>> &output) const override;
227 
228  // for documentation, see the Mapping base class
229  virtual void
230  transform(const ArrayView<const Tensor<3, dim>> & input,
231  const MappingKind kind,
232  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
233  const ArrayView<Tensor<3, spacedim>> &output) const override;
234 
244 public:
256  class InternalData : public Mapping<dim, spacedim>::InternalDataBase
257  {
258  public:
263  InternalData(const unsigned int polynomial_degree);
264 
273  void
274  initialize(const UpdateFlags update_flags,
275  const Quadrature<dim> &quadrature,
276  const unsigned int n_original_q_points);
277 
283  void
284  initialize_face(const UpdateFlags update_flags,
285  const Quadrature<dim> &quadrature,
286  const unsigned int n_original_q_points);
287 
303  void
304  compute_shape_function_values(const std::vector<Point<dim>> &unit_points);
305 
310  const double &
311  shape(const unsigned int qpoint, const unsigned int shape_nr) const;
312 
316  double &
317  shape(const unsigned int qpoint, const unsigned int shape_nr);
318 
322  const Tensor<1, dim> &
323  derivative(const unsigned int qpoint, const unsigned int shape_nr) const;
324 
329  derivative(const unsigned int qpoint, const unsigned int shape_nr);
330 
334  const Tensor<2, dim> &
335  second_derivative(const unsigned int qpoint,
336  const unsigned int shape_nr) const;
337 
342  second_derivative(const unsigned int qpoint, const unsigned int shape_nr);
343 
347  const Tensor<3, dim> &
348  third_derivative(const unsigned int qpoint,
349  const unsigned int shape_nr) const;
350 
355  third_derivative(const unsigned int qpoint, const unsigned int shape_nr);
356 
360  const Tensor<4, dim> &
361  fourth_derivative(const unsigned int qpoint,
362  const unsigned int shape_nr) const;
363 
368  fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr);
369 
373  virtual std::size_t
374  memory_consumption() const override;
375 
381  std::vector<double> shape_values;
382 
388  std::vector<Tensor<1, dim>> shape_derivatives;
389 
396  std::vector<Tensor<2, dim>> shape_second_derivatives;
397 
404  std::vector<Tensor<3, dim>> shape_third_derivatives;
405 
412  std::vector<Tensor<4, dim>> shape_fourth_derivatives;
413 
427  std::array<std::vector<Tensor<1, dim>>,
430 
435  const unsigned int polynomial_degree;
436 
446  const unsigned int n_shape_functions;
447 
448  /*
449  * The default line support points. Is used in when the shape function
450  * values are computed.
451  *
452  * The number of quadrature points depends on the degree of this
453  * class, and it matches the number of degrees of freedom of an
454  * FE_Q<1>(this->degree).
455  */
457 
465 
471 
477 
483 
489 
495 
500 
510  mutable std::vector<DerivativeForm<1, dim, spacedim>> covariant;
511 
519  mutable std::vector<DerivativeForm<1, dim, spacedim>> contravariant;
520 
524  mutable std::vector<std::vector<Tensor<1, spacedim>>> aux;
525 
530  mutable std::vector<Point<spacedim>> mapping_support_points;
531 
537 
542  mutable std::vector<double> volume_elements;
543  };
544 
545 
546  // documentation can be found in Mapping::requires_update_flags()
547  virtual UpdateFlags
548  requires_update_flags(const UpdateFlags update_flags) const override;
549 
550  // documentation can be found in Mapping::get_data()
551  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
552  get_data(const UpdateFlags, const Quadrature<dim> &quadrature) const override;
553 
554  // documentation can be found in Mapping::get_face_data()
555  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
556  get_face_data(const UpdateFlags flags,
557  const Quadrature<dim - 1> &quadrature) const override;
558 
559  // documentation can be found in Mapping::get_subface_data()
560  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
561  get_subface_data(const UpdateFlags flags,
562  const Quadrature<dim - 1> &quadrature) const override;
563 
564  // documentation can be found in Mapping::fill_fe_values()
567  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
568  const CellSimilarity::Similarity cell_similarity,
569  const Quadrature<dim> & quadrature,
570  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
572  &output_data) const override;
573 
574  // documentation can be found in Mapping::fill_fe_face_values()
575  virtual void
577  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
578  const unsigned int face_no,
579  const Quadrature<dim - 1> & quadrature,
580  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
582  &output_data) const override;
583 
584  // documentation can be found in Mapping::fill_fe_subface_values()
585  virtual void
587  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
588  const unsigned int face_no,
589  const unsigned int subface_no,
590  const Quadrature<dim - 1> & quadrature,
591  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
593  &output_data) const override;
594 
599 protected:
604  const unsigned int polynomial_degree;
605 
606  /*
607  * The default line support points. These are used when computing the
608  * location in real space of the support points on lines and quads, which
609  * are needed by the Manifold<dim,spacedim> class.
610  *
611  * The number of points depends on the degree of this class, and it matches
612  * the number of degrees of freedom of an FE_Q<1>(this->degree).
613  */
614  const std::vector<Point<1>> line_support_points;
615 
616  /*
617  * The one-dimensional polynomials defined as Lagrange polynomials from the
618  * line support points. These are used for point evaluations and match the
619  * polynomial space of an FE_Q<1>(this->degree).
620  */
621  const std::vector<Polynomials::Polynomial<double>> polynomials_1d;
622 
623  /*
624  * The numbering from the lexicographic to the hierarchical ordering used
625  * when expanding the tensor product with the mapping support points (which
626  * come in hierarchical numbers).
627  */
628  const std::vector<unsigned int> renumber_lexicographic_to_hierarchic;
629 
630  /*
631  * The support points in reference coordinates. These are used for
632  * constructing approximations of the output of
633  * compute_mapping_support_points() when evaluating the mapping on the fly,
634  * rather than going through the FEValues interface provided by
635  * InternalData.
636  *
637  * The number of points depends on the degree of this class, and it matches
638  * the number of degrees of freedom of an FE_Q<dim>(this->degree).
639  */
640  const std::vector<Point<dim>> unit_cell_support_points;
641 
661  const std::vector<Table<2, double>>
663 
677 
707  virtual std::vector<Point<spacedim>>
709  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const;
710 
715  Point<dim>
717  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
718  const Point<spacedim> & p,
719  const Point<dim> &initial_p_unit) const;
720 
734  virtual void
736  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
737  std::vector<Point<spacedim>> & a) const;
738 
753  virtual void
755  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
756  std::vector<Point<spacedim>> & a) const;
757 
758  // Make MappingQ a friend since it needs to call the fill_fe_values()
759  // functions on its MappingQGeneric(1) sub-object.
760  template <int, int>
761  friend class MappingQ;
762 
763  // Make MappingQCache a friend since it needs to call the
764  // compute_mapping_support_points() function.
765  template <int, int>
766  friend class MappingQCache;
767 };
768 
769 
770 
773 /*----------------------------------------------------------------------*/
774 
775 #ifndef DOXYGEN
776 
777 template <int dim, int spacedim>
778 inline const double &
780  const unsigned int qpoint,
781  const unsigned int shape_nr) const
782 {
783  AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
784  return shape_values[qpoint * n_shape_functions + shape_nr];
785 }
786 
787 
788 
789 template <int dim, int spacedim>
790 inline double &
792  const unsigned int shape_nr)
793 {
794  AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
795  return shape_values[qpoint * n_shape_functions + shape_nr];
796 }
797 
798 
799 template <int dim, int spacedim>
800 inline const Tensor<1, dim> &
802  const unsigned int qpoint,
803  const unsigned int shape_nr) const
804 {
805  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
806  shape_derivatives.size());
807  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
808 }
809 
810 
811 
812 template <int dim, int spacedim>
813 inline Tensor<1, dim> &
815  const unsigned int qpoint,
816  const unsigned int shape_nr)
817 {
818  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
819  shape_derivatives.size());
820  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
821 }
822 
823 
824 template <int dim, int spacedim>
825 inline const Tensor<2, dim> &
827  const unsigned int qpoint,
828  const unsigned int shape_nr) const
829 {
830  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
831  shape_second_derivatives.size());
832  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
833 }
834 
835 
836 template <int dim, int spacedim>
837 inline Tensor<2, dim> &
839  const unsigned int qpoint,
840  const unsigned int shape_nr)
841 {
842  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
843  shape_second_derivatives.size());
844  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
845 }
846 
847 template <int dim, int spacedim>
848 inline const Tensor<3, dim> &
850  const unsigned int qpoint,
851  const unsigned int shape_nr) const
852 {
853  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
854  shape_third_derivatives.size());
855  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
856 }
857 
858 
859 template <int dim, int spacedim>
860 inline Tensor<3, dim> &
862  const unsigned int qpoint,
863  const unsigned int shape_nr)
864 {
865  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
866  shape_third_derivatives.size());
867  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
868 }
869 
870 
871 template <int dim, int spacedim>
872 inline const Tensor<4, dim> &
874  const unsigned int qpoint,
875  const unsigned int shape_nr) const
876 {
877  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
878  shape_fourth_derivatives.size());
879  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
880 }
881 
882 
883 template <int dim, int spacedim>
884 inline Tensor<4, dim> &
886  const unsigned int qpoint,
887  const unsigned int shape_nr)
888 {
889  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
890  shape_fourth_derivatives.size());
891  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
892 }
893 
894 
895 
896 template <int dim, int spacedim>
897 inline bool
899 {
900  return true;
901 }
902 
903 #endif // DOXYGEN
904 
905 /* -------------- declaration of explicit specializations ------------- */
906 
907 
909 
910 #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:1691
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)
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
const std::vector< Polynomials::Polynomial< double > > polynomials_1d
QGaussLobatto< 1 > line_support_points
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
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:301
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:372
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
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
const std::vector< Point< 1 > > line_support_points
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:371
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
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