Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
mapping_q_generic.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2000 - 2019 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 
22 #include <deal.II/base/derivative_form.h>
23 #include <deal.II/base/quadrature_lib.h>
24 #include <deal.II/base/table.h>
25 #include <deal.II/base/vectorization.h>
26 
27 #include <deal.II/fe/fe_q.h>
28 #include <deal.II/fe/mapping.h>
29 
30 #include <deal.II/grid/tria_iterator.h>
31 
32 #include <deal.II/matrix_free/shape_info.h>
33 
34 #include <array>
35 #include <cmath>
36 
37 DEAL_II_NAMESPACE_OPEN
38 
39 template <int, int>
40 class MappingQ;
41 
42 
45 
46 
132 template <int dim, int spacedim = dim>
133 class MappingQGeneric : public Mapping<dim, spacedim>
134 {
135 public:
141  MappingQGeneric(const unsigned int polynomial_degree);
142 
147 
148  // for documentation, see the Mapping base class
149  virtual std::unique_ptr<Mapping<dim, spacedim>>
150  clone() const override;
151 
156  unsigned int
157  get_degree() const;
158 
163  virtual bool
164  preserves_vertex_locations() const override;
165 
171  // for documentation, see the Mapping base class
172  virtual Point<spacedim>
174  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
175  const Point<dim> &p) const override;
176 
177  // for documentation, see the Mapping base class
178  virtual Point<dim>
180  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
181  const Point<spacedim> &p) const override;
182 
192  // for documentation, see the Mapping base class
193  virtual void
194  transform(const ArrayView<const Tensor<1, dim>> & input,
195  const MappingType type,
196  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
197  const ArrayView<Tensor<1, spacedim>> &output) const override;
198 
199  // for documentation, see the Mapping base class
200  virtual void
202  const MappingType type,
203  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
204  const ArrayView<Tensor<2, spacedim>> &output) const override;
205 
206  // for documentation, see the Mapping base class
207  virtual void
208  transform(const ArrayView<const Tensor<2, dim>> & input,
209  const MappingType type,
210  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
211  const ArrayView<Tensor<2, spacedim>> &output) const override;
212 
213  // for documentation, see the Mapping base class
214  virtual void
216  const MappingType type,
217  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
218  const ArrayView<Tensor<3, spacedim>> &output) const override;
219 
220  // for documentation, see the Mapping base class
221  virtual void
222  transform(const ArrayView<const Tensor<3, dim>> & input,
223  const MappingType type,
224  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
225  const ArrayView<Tensor<3, spacedim>> &output) const override;
226 
236 public:
248  class InternalData : public Mapping<dim, spacedim>::InternalDataBase
249  {
250  public:
255  InternalData(const unsigned int polynomial_degree);
256 
265  void
266  initialize(const UpdateFlags update_flags,
267  const Quadrature<dim> &quadrature,
268  const unsigned int n_original_q_points);
269 
275  void
276  initialize_face(const UpdateFlags update_flags,
277  const Quadrature<dim> &quadrature,
278  const unsigned int n_original_q_points);
279 
295  void
296  compute_shape_function_values(const std::vector<Point<dim>> &unit_points);
297 
298 
303  const double &
304  shape(const unsigned int qpoint, const unsigned int shape_nr) const;
305 
309  double &
310  shape(const unsigned int qpoint, const unsigned int shape_nr);
311 
315  const Tensor<1, dim> &
316  derivative(const unsigned int qpoint, const unsigned int shape_nr) const;
317 
322  derivative(const unsigned int qpoint, const unsigned int shape_nr);
323 
327  const Tensor<2, dim> &
328  second_derivative(const unsigned int qpoint,
329  const unsigned int shape_nr) const;
330 
335  second_derivative(const unsigned int qpoint, const unsigned int shape_nr);
336 
340  const Tensor<3, dim> &
341  third_derivative(const unsigned int qpoint,
342  const unsigned int shape_nr) const;
343 
348  third_derivative(const unsigned int qpoint, const unsigned int shape_nr);
349 
353  const Tensor<4, dim> &
354  fourth_derivative(const unsigned int qpoint,
355  const unsigned int shape_nr) const;
356 
361  fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr);
362 
366  virtual std::size_t
367  memory_consumption() const override;
368 
374  std::vector<double> shape_values;
375 
381  std::vector<Tensor<1, dim>> shape_derivatives;
382 
389  std::vector<Tensor<2, dim>> shape_second_derivatives;
390 
397  std::vector<Tensor<3, dim>> shape_third_derivatives;
398 
405  std::vector<Tensor<4, dim>> shape_fourth_derivatives;
406 
420  std::array<std::vector<Tensor<1, dim>>,
423 
428  const unsigned int polynomial_degree;
429 
439  const unsigned int n_shape_functions;
440 
441  /*
442  * The default line support points. Is used in when the shape function
443  * values are computed.
444  *
445  * The number of quadrature points depends on the degree of this
446  * class, and it matches the number of degrees of freedom of an
447  * FE_Q<1>(this->degree).
448  */
449  QGaussLobatto<1> line_support_points;
450 
458 
464 
470 
476 
482 
488 
493 
503  mutable std::vector<DerivativeForm<1, dim, spacedim>> covariant;
504 
512  mutable std::vector<DerivativeForm<1, dim, spacedim>> contravariant;
513 
517  mutable std::vector<std::vector<Tensor<1, spacedim>>> aux;
518 
523  mutable std::vector<Point<spacedim>> mapping_support_points;
524 
530 
535  mutable std::vector<double> volume_elements;
536  };
537 
538 
539  // documentation can be found in Mapping::requires_update_flags()
540  virtual UpdateFlags
541  requires_update_flags(const UpdateFlags update_flags) const override;
542 
543  // documentation can be found in Mapping::get_data()
544  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
545  get_data(const UpdateFlags, const Quadrature<dim> &quadrature) const override;
546 
547  // documentation can be found in Mapping::get_face_data()
548  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
549  get_face_data(const UpdateFlags flags,
550  const Quadrature<dim - 1> &quadrature) const override;
551 
552  // documentation can be found in Mapping::get_subface_data()
553  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
554  get_subface_data(const UpdateFlags flags,
555  const Quadrature<dim - 1> &quadrature) const override;
556 
557  // documentation can be found in Mapping::fill_fe_values()
560  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
561  const CellSimilarity::Similarity cell_similarity,
562  const Quadrature<dim> & quadrature,
563  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
565  &output_data) const override;
566 
567  // documentation can be found in Mapping::fill_fe_face_values()
568  virtual void
570  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
571  const unsigned int face_no,
572  const Quadrature<dim - 1> & quadrature,
573  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
575  &output_data) const override;
576 
577  // documentation can be found in Mapping::fill_fe_subface_values()
578  virtual void
580  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
581  const unsigned int face_no,
582  const unsigned int subface_no,
583  const Quadrature<dim - 1> & quadrature,
584  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
586  &output_data) const override;
587 
592 protected:
597  const unsigned int polynomial_degree;
598 
599  /*
600  * The default line support points. These are used when computing
601  * the location in real space of the support points on lines and
602  * quads, which are asked to the Manifold<dim,spacedim> class.
603  *
604  * The number of quadrature points depends on the degree of this
605  * class, and it matches the number of degrees of freedom of an
606  * FE_Q<1>(this->degree).
607  */
608  QGaussLobatto<1> line_support_points;
609 
615  const std::unique_ptr<FE_Q<dim>> fe_q;
616 
636  std::vector<Table<2, double>> support_point_weights_perimeter_to_interior;
637 
651 
681  virtual std::vector<Point<spacedim>>
683  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const;
684 
689  Point<dim>
691  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
692  const Point<spacedim> & p,
693  const Point<dim> &initial_p_unit) const;
694 
708  virtual void
710  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
711  std::vector<Point<spacedim>> & a) const;
712 
727  virtual void
729  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
730  std::vector<Point<spacedim>> & a) const;
731 
736  template <int, int>
737  friend class MappingQ;
738 };
739 
740 
741 
744 /*----------------------------------------------------------------------*/
745 
746 #ifndef DOXYGEN
747 
748 template <int dim, int spacedim>
749 inline const double &
751  const unsigned int qpoint,
752  const unsigned int shape_nr) const
753 {
754  Assert(qpoint * n_shape_functions + shape_nr < shape_values.size(),
755  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
756  0,
757  shape_values.size()));
758  return shape_values[qpoint * n_shape_functions + shape_nr];
759 }
760 
761 
762 
763 template <int dim, int spacedim>
764 inline double &
766  const unsigned int shape_nr)
767 {
768  Assert(qpoint * n_shape_functions + shape_nr < shape_values.size(),
769  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
770  0,
771  shape_values.size()));
772  return shape_values[qpoint * n_shape_functions + shape_nr];
773 }
774 
775 
776 template <int dim, int spacedim>
777 inline const Tensor<1, dim> &
779  const unsigned int qpoint,
780  const unsigned int shape_nr) const
781 {
782  Assert(qpoint * n_shape_functions + shape_nr < shape_derivatives.size(),
783  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
784  0,
785  shape_derivatives.size()));
786  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
787 }
788 
789 
790 
791 template <int dim, int spacedim>
792 inline Tensor<1, dim> &
794  const unsigned int qpoint,
795  const unsigned int shape_nr)
796 {
797  Assert(qpoint * n_shape_functions + shape_nr < shape_derivatives.size(),
798  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
799  0,
800  shape_derivatives.size()));
801  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
802 }
803 
804 
805 template <int dim, int spacedim>
806 inline const Tensor<2, dim> &
808  const unsigned int qpoint,
809  const unsigned int shape_nr) const
810 {
811  Assert(qpoint * n_shape_functions + shape_nr <
812  shape_second_derivatives.size(),
813  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
814  0,
815  shape_second_derivatives.size()));
816  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
817 }
818 
819 
820 template <int dim, int spacedim>
821 inline Tensor<2, dim> &
823  const unsigned int qpoint,
824  const unsigned int shape_nr)
825 {
826  Assert(qpoint * n_shape_functions + shape_nr <
827  shape_second_derivatives.size(),
828  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
829  0,
830  shape_second_derivatives.size()));
831  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
832 }
833 
834 template <int dim, int spacedim>
835 inline const Tensor<3, dim> &
837  const unsigned int qpoint,
838  const unsigned int shape_nr) const
839 {
840  Assert(qpoint * n_shape_functions + shape_nr < shape_third_derivatives.size(),
841  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
842  0,
843  shape_third_derivatives.size()));
844  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
845 }
846 
847 
848 template <int dim, int spacedim>
849 inline Tensor<3, dim> &
851  const unsigned int qpoint,
852  const unsigned int shape_nr)
853 {
854  Assert(qpoint * n_shape_functions + shape_nr < shape_third_derivatives.size(),
855  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
856  0,
857  shape_third_derivatives.size()));
858  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
859 }
860 
861 
862 template <int dim, int spacedim>
863 inline const Tensor<4, dim> &
865  const unsigned int qpoint,
866  const unsigned int shape_nr) const
867 {
868  Assert(qpoint * n_shape_functions + shape_nr <
869  shape_fourth_derivatives.size(),
870  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
871  0,
872  shape_fourth_derivatives.size()));
873  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
874 }
875 
876 
877 template <int dim, int spacedim>
878 inline Tensor<4, dim> &
880  const unsigned int qpoint,
881  const unsigned int shape_nr)
882 {
883  Assert(qpoint * n_shape_functions + shape_nr <
884  shape_fourth_derivatives.size(),
885  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
886  0,
887  shape_fourth_derivatives.size()));
888  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
889 }
890 
891 
892 
893 template <int dim, int spacedim>
894 inline bool
896 {
897  return true;
898 }
899 
900 #endif // DOXYGEN
901 
902 /* -------------- declaration of explicit specializations ------------- */
903 
904 
905 DEAL_II_NAMESPACE_CLOSE
906 
907 #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
MappingType
Definition: mapping.h:61
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
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
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
const std::unique_ptr< FE_Q< dim > > fe_q
InternalData(const unsigned int polynomial_degree)
const unsigned int polynomial_degree
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
#define Assert(cond, exc)
Definition: exceptions.h:1407
UpdateFlags
Abstract base class for mapping classes.
Definition: dof_tools.h:57
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
virtual void transform(const ArrayView< const Tensor< 1, dim >> &input, const MappingType type, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const ArrayView< Tensor< 1, spacedim >> &output) 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
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