Reference documentation for deal.II version GIT b8135fa6eb 2023-03-25 15:55:02+00:00
\(\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.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2000 - 2022 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_h
17 #define dealii_mapping_q_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 #ifndef DOXYGEN
41 template <int, int>
42 class MappingQCache;
43 #endif
44 
110 template <int dim, int spacedim = dim>
111 class MappingQ : public Mapping<dim, spacedim>
112 {
113 public:
119  MappingQ(const unsigned int polynomial_degree);
120 
127  MappingQ(const unsigned int polynomial_degree,
128  const bool use_mapping_q_on_all_cells);
129 
133  MappingQ(const MappingQ<dim, spacedim> &mapping);
134 
135  // for documentation, see the Mapping base class
136  virtual std::unique_ptr<Mapping<dim, spacedim>>
137  clone() const override;
138 
143  unsigned int
144  get_degree() const;
145 
150  virtual bool
151  preserves_vertex_locations() const override;
152 
153  // for documentation, see the Mapping base class
154  virtual BoundingBox<spacedim>
156  &cell) const override;
157 
158  virtual bool
159  is_compatible_with(const ReferenceCell &reference_cell) const override;
160 
166  // for documentation, see the Mapping base class
167  virtual Point<spacedim>
169  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
170  const Point<dim> &p) const override;
171 
172  // for documentation, see the Mapping base class
173  virtual Point<dim>
175  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
176  const Point<spacedim> &p) const override;
177 
178  // for documentation, see the Mapping base class
179  virtual void
181  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
182  const ArrayView<const Point<spacedim>> & real_points,
183  const ArrayView<Point<dim>> &unit_points) 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 
250  void
252  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
253  const ArrayView<const Point<dim>> & unit_points,
254  const UpdateFlags update_flags,
256  &output_data) const;
257 
274  class InternalData : public Mapping<dim, spacedim>::InternalDataBase
275  {
276  public:
281  InternalData(const unsigned int polynomial_degree);
282 
291  void
292  initialize(const UpdateFlags update_flags,
293  const Quadrature<dim> &quadrature,
294  const unsigned int n_original_q_points);
295 
301  void
302  initialize_face(const UpdateFlags update_flags,
303  const Quadrature<dim> &quadrature,
304  const unsigned int n_original_q_points);
305 
321  void
322  compute_shape_function_values(const std::vector<Point<dim>> &unit_points);
323 
328  const double &
329  shape(const unsigned int qpoint, const unsigned int shape_nr) const;
330 
334  double &
335  shape(const unsigned int qpoint, const unsigned int shape_nr);
336 
340  const Tensor<1, dim> &
341  derivative(const unsigned int qpoint, const unsigned int shape_nr) const;
342 
347  derivative(const unsigned int qpoint, const unsigned int shape_nr);
348 
352  const Tensor<2, dim> &
353  second_derivative(const unsigned int qpoint,
354  const unsigned int shape_nr) const;
355 
360  second_derivative(const unsigned int qpoint, const unsigned int shape_nr);
361 
365  const Tensor<3, dim> &
366  third_derivative(const unsigned int qpoint,
367  const unsigned int shape_nr) const;
368 
373  third_derivative(const unsigned int qpoint, const unsigned int shape_nr);
374 
378  const Tensor<4, dim> &
379  fourth_derivative(const unsigned int qpoint,
380  const unsigned int shape_nr) const;
381 
386  fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr);
387 
391  virtual std::size_t
392  memory_consumption() const override;
393 
400 
407 
415 
423 
431 
445  std::array<std::vector<Tensor<1, dim>>,
448 
453  const unsigned int polynomial_degree;
454 
464  const unsigned int n_shape_functions;
465 
466  /*
467  * The default line support points. Is used in when the shape function
468  * values are computed.
469  *
470  * The number of quadrature points depends on the degree of this
471  * class, and it matches the number of degrees of freedom of an
472  * FE_Q<1>(this->degree).
473  */
475 
492  std::min<std::size_t>(VectorizedArray<double>::size(),
493  (dim <= 2 ? 2 : 4))>;
494 
501 
507 
512 
523 
532 
536  mutable std::vector<AlignedVector<Tensor<1, spacedim>>> aux;
537 
542  mutable std::vector<Point<spacedim>> mapping_support_points;
543 
549 
555  };
556 
557 protected:
558  // documentation can be found in Mapping::requires_update_flags()
559  virtual UpdateFlags
560  requires_update_flags(const UpdateFlags update_flags) const override;
561 
562  // documentation can be found in Mapping::get_data()
563  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
564  get_data(const UpdateFlags, const Quadrature<dim> &quadrature) const override;
565 
567 
568  // documentation can be found in Mapping::get_face_data()
569  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
570  get_face_data(const UpdateFlags flags,
571  const hp::QCollection<dim - 1> &quadrature) const override;
572 
573  // documentation can be found in Mapping::get_subface_data()
574  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
575  get_subface_data(const UpdateFlags flags,
576  const Quadrature<dim - 1> &quadrature) const override;
577 
578  // documentation can be found in Mapping::fill_fe_values()
581  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
582  const CellSimilarity::Similarity cell_similarity,
583  const Quadrature<dim> & quadrature,
584  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
586  &output_data) const override;
587 
589 
590  // documentation can be found in Mapping::fill_fe_face_values()
591  virtual void
593  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
594  const unsigned int face_no,
595  const hp::QCollection<dim - 1> & quadrature,
596  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
598  &output_data) const override;
599 
600  // documentation can be found in Mapping::fill_fe_subface_values()
601  virtual void
603  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
604  const unsigned int face_no,
605  const unsigned int subface_no,
606  const Quadrature<dim - 1> & quadrature,
607  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
609  &output_data) const override;
610 
611  // documentation can be found in Mapping::fill_fe_immersed_surface_values()
612  virtual void
614  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
616  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
618  &output_data) const override;
619 
628  const unsigned int polynomial_degree;
629 
630  /*
631  * The default line support points. These are used when computing the
632  * location in real space of the support points on lines and quads, which
633  * are needed by the Manifold<dim,spacedim> class.
634  *
635  * The number of points depends on the degree of this class, and it matches
636  * the number of degrees of freedom of an FE_Q<1>(this->degree).
637  */
638  const std::vector<Point<1>> line_support_points;
639 
640  /*
641  * The one-dimensional polynomials defined as Lagrange polynomials from the
642  * line support points. These are used for point evaluations and match the
643  * polynomial space of an FE_Q<1>(this->degree).
644  */
645  const std::vector<Polynomials::Polynomial<double>> polynomials_1d;
646 
647  /*
648  * The numbering from the lexicographic to the hierarchical ordering used
649  * when expanding the tensor product with the mapping support points (which
650  * come in hierarchical numbers).
651  */
652  const std::vector<unsigned int> renumber_lexicographic_to_hierarchic;
653 
654  /*
655  * The support points in reference coordinates. These are used for
656  * constructing approximations of the output of
657  * compute_mapping_support_points() when evaluating the mapping on the fly,
658  * rather than going through the FEValues interface provided by
659  * InternalData.
660  *
661  * The number of points depends on the degree of this class, and it matches
662  * the number of degrees of freedom of an FE_Q<dim>(this->degree).
663  */
664  const std::vector<Point<dim>> unit_cell_support_points;
665 
685  const std::vector<Table<2, double>>
687 
701 
731  virtual std::vector<Point<spacedim>>
733  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const;
734 
739  Point<dim>
741  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
742  const Point<spacedim> & p,
743  const Point<dim> &initial_p_unit) const;
744 
758  virtual void
760  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
761  std::vector<Point<spacedim>> & a) const;
762 
777  virtual void
779  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
780  std::vector<Point<spacedim>> & a) const;
781 
782  // Make MappingQCache a friend since it needs to call the
783  // compute_mapping_support_points() function.
784  template <int, int>
785  friend class MappingQCache;
786 };
787 
788 
789 
795 template <int dim, int spacedim = dim>
797 
801 /*----------------------------------------------------------------------*/
802 
803 #ifndef DOXYGEN
804 
805 template <int dim, int spacedim>
806 inline const double &
807 MappingQ<dim, spacedim>::InternalData::shape(const unsigned int qpoint,
808  const unsigned int shape_nr) const
809 {
810  AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
811  return shape_values[qpoint * n_shape_functions + shape_nr];
812 }
813 
814 
815 
816 template <int dim, int spacedim>
817 inline double &
818 MappingQ<dim, spacedim>::InternalData::shape(const unsigned int qpoint,
819  const unsigned int shape_nr)
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 template <int dim, int spacedim>
827 inline const Tensor<1, dim> &
829  const unsigned int qpoint,
830  const unsigned int shape_nr) const
831 {
832  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
833  shape_derivatives.size());
834  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
835 }
836 
837 
838 
839 template <int dim, int spacedim>
840 inline Tensor<1, dim> &
841 MappingQ<dim, spacedim>::InternalData::derivative(const unsigned int qpoint,
842  const unsigned int shape_nr)
843 {
844  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
845  shape_derivatives.size());
846  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
847 }
848 
849 
850 template <int dim, int spacedim>
851 inline const Tensor<2, dim> &
853  const unsigned int qpoint,
854  const unsigned int shape_nr) const
855 {
856  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
857  shape_second_derivatives.size());
858  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
859 }
860 
861 
862 template <int dim, int spacedim>
863 inline Tensor<2, dim> &
865  const unsigned int qpoint,
866  const unsigned int shape_nr)
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 template <int dim, int spacedim>
874 inline const Tensor<3, dim> &
876  const unsigned int qpoint,
877  const unsigned int shape_nr) const
878 {
879  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
880  shape_third_derivatives.size());
881  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
882 }
883 
884 
885 template <int dim, int spacedim>
886 inline Tensor<3, dim> &
888  const unsigned int qpoint,
889  const unsigned int shape_nr)
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 const Tensor<4, dim> &
900  const unsigned int qpoint,
901  const unsigned int shape_nr) const
902 {
903  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
904  shape_fourth_derivatives.size());
905  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
906 }
907 
908 
909 template <int dim, int spacedim>
910 inline Tensor<4, dim> &
912  const unsigned int qpoint,
913  const unsigned int shape_nr)
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 
922 template <int dim, int spacedim>
923 inline bool
925 {
926  return true;
927 }
928 
929 #endif // DOXYGEN
930 
931 /* -------------- declaration of explicit specializations ------------- */
932 
933 
935 
936 #endif
size_type size() const
Triangulation< dim, spacedim >::cell_iterator cell_of_current_support_points
Definition: mapping_q.h:548
AlignedVector< DerivativeForm< 1, dim, spacedim > > covariant
Definition: mapping_q.h:522
QGaussLobatto< 1 > line_support_points
Definition: mapping_q.h:474
AlignedVector< Tensor< 3, dim > > shape_third_derivatives
Definition: mapping_q.h:422
Tensor< 2, dim > & second_derivative(const unsigned int qpoint, const unsigned int shape_nr)
virtual std::size_t memory_consumption() const override
Definition: mapping_q.cc:63
std::vector< Point< spacedim > > mapping_support_points
Definition: mapping_q.h:542
const Tensor< 1, dim > & derivative(const unsigned int qpoint, const unsigned int shape_nr) const
AlignedVector< Tensor< 4, dim > > shape_fourth_derivatives
Definition: mapping_q.h:430
const double & shape(const unsigned int qpoint, const unsigned int shape_nr) const
const unsigned int polynomial_degree
Definition: mapping_q.h:453
AlignedVector< double > shape_values
Definition: mapping_q.h:399
std::array< std::vector< Tensor< 1, dim > >, GeometryInfo< dim >::faces_per_cell *(dim - 1)> unit_tangentials
Definition: mapping_q.h:447
AlignedVector< DerivativeForm< 1, dim, spacedim > > contravariant
Definition: mapping_q.h:531
AlignedVector< Tensor< 1, dim > > shape_derivatives
Definition: mapping_q.h:406
const Tensor< 3, dim > & third_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
const Tensor< 2, dim > & second_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
AlignedVector< Tensor< 2, dim > > shape_second_derivatives
Definition: mapping_q.h:414
std::vector< AlignedVector< Tensor< 1, spacedim > > > aux
Definition: mapping_q.h:536
double & shape(const unsigned int qpoint, const unsigned int shape_nr)
void initialize_face(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
Definition: mapping_q.cc:214
InternalData(const unsigned int polynomial_degree)
Definition: mapping_q.cc:51
void compute_shape_function_values(const std::vector< Point< dim >> &unit_points)
Definition: mapping_q.cc:270
AlignedVector< double > volume_elements
Definition: mapping_q.h:554
void initialize(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
Definition: mapping_q.cc:84
Tensor< 4, dim > & fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr)
Tensor< 3, dim > & third_derivative(const unsigned int qpoint, const unsigned int shape_nr)
Tensor< 1, dim > & derivative(const unsigned int qpoint, const unsigned int shape_nr)
const Tensor< 4, dim > & fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
const unsigned int n_shape_functions
Definition: mapping_q.h:464
internal::MatrixFreeFunctions::ShapeInfo< VectorizedArrayType > shape_info
Definition: mapping_q.h:500
AlignedVector< VectorizedArrayType > scratch
Definition: mapping_q.h:506
virtual void add_line_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim >> &a) const
Definition: mapping_q.cc:1637
const std::vector< unsigned int > renumber_lexicographic_to_hierarchic
Definition: mapping_q.h:652
const Table< 2, double > support_point_weights_cell
Definition: mapping_q.h:700
virtual std::vector< Point< spacedim > > compute_mapping_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const
Definition: mapping_q.cc:1814
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const override
Definition: mapping_q.cc:889
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
Definition: mapping_q.cc:1398
virtual BoundingBox< spacedim > get_bounding_box(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const override
Definition: mapping_q.cc:1913
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags flags, const Quadrature< dim - 1 > &quadrature) const override
Definition: mapping_q.cc:925
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
Definition: mapping_q.cc:1215
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
Definition: mapping_q.cc:944
virtual void fill_fe_immersed_surface_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const NonMatching::ImmersedSurfaceQuadrature< dim > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
Definition: mapping_q.cc:1265
const unsigned int polynomial_degree
Definition: mapping_q.h:628
virtual bool preserves_vertex_locations() const override
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
Definition: mapping_q.cc:1499
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_face_data(const UpdateFlags flags, const hp::QCollection< dim - 1 > &quadrature) const override
Definition: mapping_q.cc:904
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
Definition: mapping_q.cc:487
const std::vector< Table< 2, double > > support_point_weights_perimeter_to_interior
Definition: mapping_q.h:686
const std::vector< Point< 1 > > line_support_points
Definition: mapping_q.h:638
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const override
Definition: mapping_q.cc:452
const std::vector< Polynomials::Polynomial< double > > polynomials_1d
Definition: mapping_q.h:645
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
Definition: mapping_q.cc:434
const std::vector< Point< dim > > unit_cell_support_points
Definition: mapping_q.h:664
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
Definition: mapping_q.cc:751
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const override
Definition: mapping_q.cc:630
MappingQ(const unsigned int polynomial_degree)
Definition: mapping_q.cc:361
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
Definition: mapping_q.cc:833
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
Definition: mapping_q.cc:1166
virtual bool is_compatible_with(const ReferenceCell &reference_cell) const override
Definition: mapping_q.cc:1923
virtual void add_quad_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim >> &a) const
Definition: mapping_q.cc:1803
unsigned int get_degree() const
Definition: mapping_q.cc:443
Abstract base class for mapping classes.
Definition: mapping.h:314
#define DEAL_II_DEPRECATED
Definition: config.h:174
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:474
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:475
#define AssertIndexRange(index, range)
Definition: exceptions.h:1827
UpdateFlags
MappingKind
Definition: mapping.h:75
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)