Reference documentation for deal.II version GIT 969cfe8cd8 2022-06-28 20:20:01+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 
47 
48 
109 template <int dim, int spacedim = dim>
110 class MappingQ : public Mapping<dim, spacedim>
111 {
112 public:
118  MappingQ(const unsigned int polynomial_degree);
119 
126  MappingQ(const unsigned int polynomial_degree,
127  const bool use_mapping_q_on_all_cells);
128 
132  MappingQ(const MappingQ<dim, spacedim> &mapping);
133 
134  // for documentation, see the Mapping base class
135  virtual std::unique_ptr<Mapping<dim, spacedim>>
136  clone() const override;
137 
142  unsigned int
143  get_degree() const;
144 
149  virtual bool
150  preserves_vertex_locations() const override;
151 
152  // for documentation, see the Mapping base class
153  virtual BoundingBox<spacedim>
155  &cell) const override;
156 
157  virtual bool
158  is_compatible_with(const ReferenceCell &reference_cell) const override;
159 
165  // for documentation, see the Mapping base class
166  virtual Point<spacedim>
168  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
169  const Point<dim> &p) const override;
170 
171  // for documentation, see the Mapping base class
172  virtual Point<dim>
174  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
175  const Point<spacedim> &p) const override;
176 
177  // for documentation, see the Mapping base class
178  virtual void
180  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
181  const ArrayView<const Point<spacedim>> & real_points,
182  const ArrayView<Point<dim>> &unit_points) const override;
183 
193  // for documentation, see the Mapping base class
194  virtual void
195  transform(const ArrayView<const Tensor<1, dim>> & input,
196  const MappingKind kind,
197  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
198  const ArrayView<Tensor<1, spacedim>> &output) const override;
199 
200  // for documentation, see the Mapping base class
201  virtual void
203  const MappingKind kind,
204  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
205  const ArrayView<Tensor<2, spacedim>> &output) const override;
206 
207  // for documentation, see the Mapping base class
208  virtual void
209  transform(const ArrayView<const Tensor<2, dim>> & input,
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
217  const MappingKind kind,
218  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
219  const ArrayView<Tensor<3, spacedim>> &output) const override;
220 
221  // for documentation, see the Mapping base class
222  virtual void
223  transform(const ArrayView<const Tensor<3, dim>> & input,
224  const MappingKind kind,
225  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
226  const ArrayView<Tensor<3, spacedim>> &output) const override;
227 
249  void
251  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
252  const ArrayView<const Point<dim>> & unit_points,
253  const UpdateFlags update_flags,
255  &output_data) const;
256 
273  class InternalData : public Mapping<dim, spacedim>::InternalDataBase
274  {
275  public:
280  InternalData(const unsigned int polynomial_degree);
281 
290  void
291  initialize(const UpdateFlags update_flags,
292  const Quadrature<dim> &quadrature,
293  const unsigned int n_original_q_points);
294 
300  void
301  initialize_face(const UpdateFlags update_flags,
302  const Quadrature<dim> &quadrature,
303  const unsigned int n_original_q_points);
304 
320  void
321  compute_shape_function_values(const std::vector<Point<dim>> &unit_points);
322 
327  const double &
328  shape(const unsigned int qpoint, const unsigned int shape_nr) const;
329 
333  double &
334  shape(const unsigned int qpoint, const unsigned int shape_nr);
335 
339  const Tensor<1, dim> &
340  derivative(const unsigned int qpoint, const unsigned int shape_nr) const;
341 
346  derivative(const unsigned int qpoint, const unsigned int shape_nr);
347 
351  const Tensor<2, dim> &
352  second_derivative(const unsigned int qpoint,
353  const unsigned int shape_nr) const;
354 
359  second_derivative(const unsigned int qpoint, const unsigned int shape_nr);
360 
364  const Tensor<3, dim> &
365  third_derivative(const unsigned int qpoint,
366  const unsigned int shape_nr) const;
367 
372  third_derivative(const unsigned int qpoint, const unsigned int shape_nr);
373 
377  const Tensor<4, dim> &
378  fourth_derivative(const unsigned int qpoint,
379  const unsigned int shape_nr) const;
380 
385  fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr);
386 
390  virtual std::size_t
391  memory_consumption() const override;
392 
399 
406 
414 
422 
430 
444  std::array<std::vector<Tensor<1, dim>>,
447 
452  const unsigned int polynomial_degree;
453 
463  const unsigned int n_shape_functions;
464 
465  /*
466  * The default line support points. Is used in when the shape function
467  * values are computed.
468  *
469  * The number of quadrature points depends on the degree of this
470  * class, and it matches the number of degrees of freedom of an
471  * FE_Q<1>(this->degree).
472  */
474 
491  std::min<std::size_t>(VectorizedArray<double>::size(),
492  (dim <= 2 ? 2 : 4))>;
493 
500 
506 
511 
522 
531 
535  mutable std::vector<AlignedVector<Tensor<1, spacedim>>> aux;
536 
541  mutable std::vector<Point<spacedim>> mapping_support_points;
542 
548 
554  };
555 
556 protected:
557  // documentation can be found in Mapping::requires_update_flags()
558  virtual UpdateFlags
559  requires_update_flags(const UpdateFlags update_flags) const override;
560 
561  // documentation can be found in Mapping::get_data()
562  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
563  get_data(const UpdateFlags, const Quadrature<dim> &quadrature) const override;
564 
566 
567  // documentation can be found in Mapping::get_face_data()
568  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
569  get_face_data(const UpdateFlags flags,
570  const hp::QCollection<dim - 1> &quadrature) const override;
571 
572  // documentation can be found in Mapping::get_subface_data()
573  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
574  get_subface_data(const UpdateFlags flags,
575  const Quadrature<dim - 1> &quadrature) const override;
576 
577  // documentation can be found in Mapping::fill_fe_values()
580  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
581  const CellSimilarity::Similarity cell_similarity,
582  const Quadrature<dim> & quadrature,
583  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
585  &output_data) const override;
586 
588 
589  // documentation can be found in Mapping::fill_fe_face_values()
590  virtual void
592  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
593  const unsigned int face_no,
594  const hp::QCollection<dim - 1> & quadrature,
595  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
597  &output_data) const override;
598 
599  // documentation can be found in Mapping::fill_fe_subface_values()
600  virtual void
602  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
603  const unsigned int face_no,
604  const unsigned int subface_no,
605  const Quadrature<dim - 1> & quadrature,
606  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
608  &output_data) const override;
609 
610  // documentation can be found in Mapping::fill_fe_immersed_surface_values()
611  virtual void
613  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
615  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
617  &output_data) const override;
618 
627  const unsigned int polynomial_degree;
628 
629  /*
630  * The default line support points. These are used when computing the
631  * location in real space of the support points on lines and quads, which
632  * are needed by the Manifold<dim,spacedim> class.
633  *
634  * The number of points depends on the degree of this class, and it matches
635  * the number of degrees of freedom of an FE_Q<1>(this->degree).
636  */
637  const std::vector<Point<1>> line_support_points;
638 
639  /*
640  * The one-dimensional polynomials defined as Lagrange polynomials from the
641  * line support points. These are used for point evaluations and match the
642  * polynomial space of an FE_Q<1>(this->degree).
643  */
644  const std::vector<Polynomials::Polynomial<double>> polynomials_1d;
645 
646  /*
647  * The numbering from the lexicographic to the hierarchical ordering used
648  * when expanding the tensor product with the mapping support points (which
649  * come in hierarchical numbers).
650  */
651  const std::vector<unsigned int> renumber_lexicographic_to_hierarchic;
652 
653  /*
654  * The support points in reference coordinates. These are used for
655  * constructing approximations of the output of
656  * compute_mapping_support_points() when evaluating the mapping on the fly,
657  * rather than going through the FEValues interface provided by
658  * InternalData.
659  *
660  * The number of points depends on the degree of this class, and it matches
661  * the number of degrees of freedom of an FE_Q<dim>(this->degree).
662  */
663  const std::vector<Point<dim>> unit_cell_support_points;
664 
684  const std::vector<Table<2, double>>
686 
700 
730  virtual std::vector<Point<spacedim>>
732  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const;
733 
738  Point<dim>
740  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
741  const Point<spacedim> & p,
742  const Point<dim> &initial_p_unit) const;
743 
757  virtual void
759  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
760  std::vector<Point<spacedim>> & a) const;
761 
776  virtual void
778  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
779  std::vector<Point<spacedim>> & a) const;
780 
781  // Make MappingQCache a friend since it needs to call the
782  // compute_mapping_support_points() function.
783  template <int, int>
784  friend class MappingQCache;
785 };
786 
787 
788 
794 template <int dim, int spacedim = dim>
796 
800 /*----------------------------------------------------------------------*/
801 
802 #ifndef DOXYGEN
803 
804 template <int dim, int spacedim>
805 inline const double &
806 MappingQ<dim, spacedim>::InternalData::shape(const unsigned int qpoint,
807  const unsigned int shape_nr) const
808 {
809  AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
810  return shape_values[qpoint * n_shape_functions + shape_nr];
811 }
812 
813 
814 
815 template <int dim, int spacedim>
816 inline double &
817 MappingQ<dim, spacedim>::InternalData::shape(const unsigned int qpoint,
818  const unsigned int shape_nr)
819 {
820  AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
821  return shape_values[qpoint * n_shape_functions + shape_nr];
822 }
823 
824 
825 template <int dim, int spacedim>
826 inline const Tensor<1, dim> &
828  const unsigned int qpoint,
829  const unsigned int shape_nr) const
830 {
831  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
832  shape_derivatives.size());
833  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
834 }
835 
836 
837 
838 template <int dim, int spacedim>
839 inline Tensor<1, dim> &
840 MappingQ<dim, spacedim>::InternalData::derivative(const unsigned int qpoint,
841  const unsigned int shape_nr)
842 {
843  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
844  shape_derivatives.size());
845  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
846 }
847 
848 
849 template <int dim, int spacedim>
850 inline const Tensor<2, dim> &
852  const unsigned int qpoint,
853  const unsigned int shape_nr) const
854 {
855  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
856  shape_second_derivatives.size());
857  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
858 }
859 
860 
861 template <int dim, int spacedim>
862 inline Tensor<2, dim> &
864  const unsigned int qpoint,
865  const unsigned int shape_nr)
866 {
867  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
868  shape_second_derivatives.size());
869  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
870 }
871 
872 template <int dim, int spacedim>
873 inline const Tensor<3, dim> &
875  const unsigned int qpoint,
876  const unsigned int shape_nr) const
877 {
878  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
879  shape_third_derivatives.size());
880  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
881 }
882 
883 
884 template <int dim, int spacedim>
885 inline Tensor<3, dim> &
887  const unsigned int qpoint,
888  const unsigned int shape_nr)
889 {
890  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
891  shape_third_derivatives.size());
892  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
893 }
894 
895 
896 template <int dim, int spacedim>
897 inline const Tensor<4, dim> &
899  const unsigned int qpoint,
900  const unsigned int shape_nr) const
901 {
902  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
903  shape_fourth_derivatives.size());
904  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
905 }
906 
907 
908 template <int dim, int spacedim>
909 inline Tensor<4, dim> &
911  const unsigned int qpoint,
912  const unsigned int shape_nr)
913 {
914  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
915  shape_fourth_derivatives.size());
916  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
917 }
918 
919 
920 
921 template <int dim, int spacedim>
922 inline bool
924 {
925  return true;
926 }
927 
928 #endif // DOXYGEN
929 
930 /* -------------- declaration of explicit specializations ------------- */
931 
932 
934 
935 #endif
size_type size() const
Triangulation< dim, spacedim >::cell_iterator cell_of_current_support_points
Definition: mapping_q.h:547
AlignedVector< DerivativeForm< 1, dim, spacedim > > covariant
Definition: mapping_q.h:521
QGaussLobatto< 1 > line_support_points
Definition: mapping_q.h:473
AlignedVector< Tensor< 3, dim > > shape_third_derivatives
Definition: mapping_q.h:421
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:64
std::vector< Point< spacedim > > mapping_support_points
Definition: mapping_q.h:541
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:429
const double & shape(const unsigned int qpoint, const unsigned int shape_nr) const
const unsigned int polynomial_degree
Definition: mapping_q.h:452
AlignedVector< double > shape_values
Definition: mapping_q.h:398
std::array< std::vector< Tensor< 1, dim > >, GeometryInfo< dim >::faces_per_cell *(dim - 1)> unit_tangentials
Definition: mapping_q.h:446
AlignedVector< DerivativeForm< 1, dim, spacedim > > contravariant
Definition: mapping_q.h:530
AlignedVector< Tensor< 1, dim > > shape_derivatives
Definition: mapping_q.h:405
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:413
std::vector< AlignedVector< Tensor< 1, spacedim > > > aux
Definition: mapping_q.h:535
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:215
InternalData(const unsigned int polynomial_degree)
Definition: mapping_q.cc:52
void compute_shape_function_values(const std::vector< Point< dim >> &unit_points)
Definition: mapping_q.cc:271
AlignedVector< double > volume_elements
Definition: mapping_q.h:553
void initialize(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
Definition: mapping_q.cc:85
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:463
internal::MatrixFreeFunctions::ShapeInfo< VectorizedArrayType > shape_info
Definition: mapping_q.h:499
AlignedVector< VectorizedArrayType > scratch
Definition: mapping_q.h:505
virtual void add_line_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim >> &a) const
Definition: mapping_q.cc:1636
const std::vector< unsigned int > renumber_lexicographic_to_hierarchic
Definition: mapping_q.h:651
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:945
const Table< 2, double > support_point_weights_cell
Definition: mapping_q.h:699
virtual std::vector< Point< spacedim > > compute_mapping_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const
Definition: mapping_q.cc:1813
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const override
Definition: mapping_q.cc:890
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:1397
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:1214
virtual BoundingBox< spacedim > get_bounding_box(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const override
Definition: mapping_q.cc:1912
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:926
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:1165
const unsigned int polynomial_degree
Definition: mapping_q.h:627
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:1498
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:905
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:488
const std::vector< Table< 2, double > > support_point_weights_perimeter_to_interior
Definition: mapping_q.h:685
const std::vector< Point< 1 > > line_support_points
Definition: mapping_q.h:637
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:453
const std::vector< Polynomials::Polynomial< double > > polynomials_1d
Definition: mapping_q.h:644
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
Definition: mapping_q.cc:435
const std::vector< Point< dim > > unit_cell_support_points
Definition: mapping_q.h:663
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:752
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:631
MappingQ(const unsigned int polynomial_degree)
Definition: mapping_q.cc:362
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
Definition: mapping_q.cc:834
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:1264
virtual bool is_compatible_with(const ReferenceCell &reference_cell) const override
Definition: mapping_q.cc:1922
virtual void add_quad_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim >> &a) const
Definition: mapping_q.cc:1802
unsigned int get_degree() const
Definition: mapping_q.cc:444
Abstract base class for mapping classes.
Definition: mapping.h:311
#define DEAL_II_DEPRECATED
Definition: config.h:164
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
UpdateFlags
#define AssertIndexRange(index, range)
Definition: exceptions.h:1732
MappingKind
Definition: mapping.h:72
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)