Reference documentation for deal.II version Git 5f32ffafa8 2019-10-18 00:26:50 -0600
\(\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 template <int, int>
43 class MappingQCache;
44 
45 
48 
49 
135 template <int dim, int spacedim = dim>
136 class MappingQGeneric : public Mapping<dim, spacedim>
137 {
138 public:
144  MappingQGeneric(const unsigned int polynomial_degree);
145 
150 
151  // for documentation, see the Mapping base class
152  virtual std::unique_ptr<Mapping<dim, spacedim>>
153  clone() const override;
154 
159  unsigned int
160  get_degree() const;
161 
166  virtual bool
167  preserves_vertex_locations() const override;
168 
174  // for documentation, see the Mapping base class
175  virtual Point<spacedim>
177  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
178  const Point<dim> &p) const override;
179 
180  // for documentation, see the Mapping base class
181  virtual Point<dim>
183  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
184  const Point<spacedim> &p) 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 
239 public:
251  class InternalData : public Mapping<dim, spacedim>::InternalDataBase
252  {
253  public:
258  InternalData(const unsigned int polynomial_degree);
259 
268  void
269  initialize(const UpdateFlags update_flags,
270  const Quadrature<dim> &quadrature,
271  const unsigned int n_original_q_points);
272 
278  void
279  initialize_face(const UpdateFlags update_flags,
280  const Quadrature<dim> &quadrature,
281  const unsigned int n_original_q_points);
282 
298  void
299  compute_shape_function_values(const std::vector<Point<dim>> &unit_points);
300 
301 
306  const double &
307  shape(const unsigned int qpoint, const unsigned int shape_nr) const;
308 
312  double &
313  shape(const unsigned int qpoint, const unsigned int shape_nr);
314 
318  const Tensor<1, dim> &
319  derivative(const unsigned int qpoint, const unsigned int shape_nr) const;
320 
325  derivative(const unsigned int qpoint, const unsigned int shape_nr);
326 
330  const Tensor<2, dim> &
331  second_derivative(const unsigned int qpoint,
332  const unsigned int shape_nr) const;
333 
338  second_derivative(const unsigned int qpoint, const unsigned int shape_nr);
339 
343  const Tensor<3, dim> &
344  third_derivative(const unsigned int qpoint,
345  const unsigned int shape_nr) const;
346 
351  third_derivative(const unsigned int qpoint, const unsigned int shape_nr);
352 
356  const Tensor<4, dim> &
357  fourth_derivative(const unsigned int qpoint,
358  const unsigned int shape_nr) const;
359 
364  fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr);
365 
369  virtual std::size_t
370  memory_consumption() const override;
371 
377  std::vector<double> shape_values;
378 
384  std::vector<Tensor<1, dim>> shape_derivatives;
385 
392  std::vector<Tensor<2, dim>> shape_second_derivatives;
393 
400  std::vector<Tensor<3, dim>> shape_third_derivatives;
401 
408  std::vector<Tensor<4, dim>> shape_fourth_derivatives;
409 
423  std::array<std::vector<Tensor<1, dim>>,
426 
431  const unsigned int polynomial_degree;
432 
442  const unsigned int n_shape_functions;
443 
444  /*
445  * The default line support points. Is used in when the shape function
446  * values are computed.
447  *
448  * The number of quadrature points depends on the degree of this
449  * class, and it matches the number of degrees of freedom of an
450  * FE_Q<1>(this->degree).
451  */
452  QGaussLobatto<1> line_support_points;
453 
461 
467 
473 
479 
485 
491 
496 
506  mutable std::vector<DerivativeForm<1, dim, spacedim>> covariant;
507 
515  mutable std::vector<DerivativeForm<1, dim, spacedim>> contravariant;
516 
520  mutable std::vector<std::vector<Tensor<1, spacedim>>> aux;
521 
526  mutable std::vector<Point<spacedim>> mapping_support_points;
527 
533 
538  mutable std::vector<double> volume_elements;
539  };
540 
541 
542  // documentation can be found in Mapping::requires_update_flags()
543  virtual UpdateFlags
544  requires_update_flags(const UpdateFlags update_flags) const override;
545 
546  // documentation can be found in Mapping::get_data()
547  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
548  get_data(const UpdateFlags, const Quadrature<dim> &quadrature) const override;
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 Quadrature<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 
570  // documentation can be found in Mapping::fill_fe_face_values()
571  virtual void
573  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
574  const unsigned int face_no,
575  const Quadrature<dim - 1> & quadrature,
576  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
578  &output_data) const override;
579 
580  // documentation can be found in Mapping::fill_fe_subface_values()
581  virtual void
583  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
584  const unsigned int face_no,
585  const unsigned int subface_no,
586  const Quadrature<dim - 1> & quadrature,
587  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
589  &output_data) const override;
590 
595 protected:
600  const unsigned int polynomial_degree;
601 
602  /*
603  * The default line support points. These are used when computing
604  * the location in real space of the support points on lines and
605  * quads, which are asked to the Manifold<dim,spacedim> class.
606  *
607  * The number of quadrature points depends on the degree of this
608  * class, and it matches the number of degrees of freedom of an
609  * FE_Q<1>(this->degree).
610  */
611  QGaussLobatto<1> line_support_points;
612 
618  const std::unique_ptr<FE_Q<dim>> fe_q;
619 
639  std::vector<Table<2, double>> support_point_weights_perimeter_to_interior;
640 
654 
684  virtual std::vector<Point<spacedim>>
686  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const;
687 
692  Point<dim>
694  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
695  const Point<spacedim> & p,
696  const Point<dim> &initial_p_unit) const;
697 
711  virtual void
713  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
714  std::vector<Point<spacedim>> & a) const;
715 
730  virtual void
732  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
733  std::vector<Point<spacedim>> & a) const;
734 
735  // Make MappingQ a friend since it needs to call the fill_fe_values()
736  // functions on its MappingQGeneric(1) sub-object.
737  template <int, int>
738  friend class MappingQ;
739 
740  // Make MappingQCache a friend since it needs to call the
741  // compute_mapping_support_points() function.
742  template <int, int>
743  friend class MappingQCache;
744 };
745 
746 
747 
750 /*----------------------------------------------------------------------*/
751 
752 #ifndef DOXYGEN
753 
754 template <int dim, int spacedim>
755 inline const double &
757  const unsigned int qpoint,
758  const unsigned int shape_nr) const
759 {
760  Assert(qpoint * n_shape_functions + shape_nr < shape_values.size(),
761  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
762  0,
763  shape_values.size()));
764  return shape_values[qpoint * n_shape_functions + shape_nr];
765 }
766 
767 
768 
769 template <int dim, int spacedim>
770 inline double &
772  const unsigned int shape_nr)
773 {
774  Assert(qpoint * n_shape_functions + shape_nr < shape_values.size(),
775  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
776  0,
777  shape_values.size()));
778  return shape_values[qpoint * n_shape_functions + shape_nr];
779 }
780 
781 
782 template <int dim, int spacedim>
783 inline const Tensor<1, dim> &
785  const unsigned int qpoint,
786  const unsigned int shape_nr) const
787 {
788  Assert(qpoint * n_shape_functions + shape_nr < shape_derivatives.size(),
789  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
790  0,
791  shape_derivatives.size()));
792  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
793 }
794 
795 
796 
797 template <int dim, int spacedim>
798 inline Tensor<1, dim> &
800  const unsigned int qpoint,
801  const unsigned int shape_nr)
802 {
803  Assert(qpoint * n_shape_functions + shape_nr < shape_derivatives.size(),
804  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
805  0,
806  shape_derivatives.size()));
807  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
808 }
809 
810 
811 template <int dim, int spacedim>
812 inline const Tensor<2, dim> &
814  const unsigned int qpoint,
815  const unsigned int shape_nr) const
816 {
817  Assert(qpoint * n_shape_functions + shape_nr <
819  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
820  0,
821  shape_second_derivatives.size()));
822  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
823 }
824 
825 
826 template <int dim, int spacedim>
827 inline Tensor<2, dim> &
829  const unsigned int qpoint,
830  const unsigned int shape_nr)
831 {
832  Assert(qpoint * n_shape_functions + shape_nr <
834  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
835  0,
836  shape_second_derivatives.size()));
837  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
838 }
839 
840 template <int dim, int spacedim>
841 inline const Tensor<3, dim> &
843  const unsigned int qpoint,
844  const unsigned int shape_nr) const
845 {
846  Assert(qpoint * n_shape_functions + shape_nr < shape_third_derivatives.size(),
847  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
848  0,
849  shape_third_derivatives.size()));
850  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
851 }
852 
853 
854 template <int dim, int spacedim>
855 inline Tensor<3, dim> &
857  const unsigned int qpoint,
858  const unsigned int shape_nr)
859 {
860  Assert(qpoint * n_shape_functions + shape_nr < shape_third_derivatives.size(),
861  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
862  0,
863  shape_third_derivatives.size()));
864  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
865 }
866 
867 
868 template <int dim, int spacedim>
869 inline const Tensor<4, dim> &
871  const unsigned int qpoint,
872  const unsigned int shape_nr) const
873 {
874  Assert(qpoint * n_shape_functions + shape_nr <
876  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
877  0,
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  Assert(qpoint * n_shape_functions + shape_nr <
891  ExcIndexRange(qpoint * n_shape_functions + shape_nr,
892  0,
893  shape_fourth_derivatives.size()));
894  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
895 }
896 
897 
898 
899 template <int dim, int spacedim>
900 inline bool
902 {
903  return true;
904 }
905 
906 #endif // DOXYGEN
907 
908 /* -------------- declaration of explicit specializations ------------- */
909 
910 
911 DEAL_II_NAMESPACE_CLOSE
912 
913 #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
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)
MappingKind
Definition: mapping.h:62
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
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
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
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