Reference documentation for deal.II version GIT c9976103bc 2022-12-09 17:30: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_fe.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2020 - 2021 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_fe_h
17 #define dealii_mapping_fe_h
18 
19 
20 #include <deal.II/base/config.h>
21 
24 #include <deal.II/base/table.h>
26 
27 #include <deal.II/fe/mapping.h>
28 
30 
31 #include <array>
32 #include <cmath>
33 
35 
36 
58 template <int dim, int spacedim = dim>
59 class MappingFE : public Mapping<dim, spacedim>
60 {
61 public:
65  explicit MappingFE(const FiniteElement<dim, spacedim> &fe);
66 
70  MappingFE(const MappingFE<dim, spacedim> &mapping);
71 
72  // for documentation, see the Mapping base class
73  virtual std::unique_ptr<Mapping<dim, spacedim>>
74  clone() const override;
75 
80  unsigned int
81  get_degree() const;
82 
83  // for documentation, see the Mapping base class
84  virtual BoundingBox<spacedim>
86  &cell) const override;
87 
88 
89  virtual bool
90  is_compatible_with(const ReferenceCell &reference_cell) const override;
91 
96  virtual bool
97  preserves_vertex_locations() const override;
98 
104  // for documentation, see the Mapping base class
105  virtual Point<spacedim>
107  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
108  const Point<dim> &p) const override;
109 
110  // for documentation, see the Mapping base class
111  virtual Point<dim>
113  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
114  const Point<spacedim> &p) const override;
115 
125  // for documentation, see the Mapping base class
126  virtual void
127  transform(const ArrayView<const Tensor<1, dim>> & input,
128  const MappingKind kind,
129  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
130  const ArrayView<Tensor<1, spacedim>> &output) const override;
131 
132  // for documentation, see the Mapping base class
133  virtual void
135  const MappingKind kind,
136  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
137  const ArrayView<Tensor<2, spacedim>> &output) const override;
138 
139  // for documentation, see the Mapping base class
140  virtual void
141  transform(const ArrayView<const Tensor<2, dim>> & input,
142  const MappingKind kind,
143  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
144  const ArrayView<Tensor<2, spacedim>> &output) const override;
145 
146  // for documentation, see the Mapping base class
147  virtual void
149  const MappingKind kind,
150  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
151  const ArrayView<Tensor<3, spacedim>> &output) const override;
152 
153  // for documentation, see the Mapping base class
154  virtual void
155  transform(const ArrayView<const Tensor<3, dim>> & input,
156  const MappingKind kind,
157  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
158  const ArrayView<Tensor<3, spacedim>> &output) const override;
159 
180  class InternalData : public Mapping<dim, spacedim>::InternalDataBase
181  {
182  public:
187 
196  void
197  initialize(const UpdateFlags update_flags,
198  const Quadrature<dim> &quadrature,
199  const unsigned int n_original_q_points);
200 
206  void
207  initialize_face(const UpdateFlags update_flags,
208  const Quadrature<dim> &quadrature,
209  const unsigned int n_original_q_points);
210 
215  void
216  compute_shape_function_values(const std::vector<Point<dim>> &unit_points);
217 
218 
223  const double &
224  shape(const unsigned int qpoint, const unsigned int shape_nr) const;
225 
229  double &
230  shape(const unsigned int qpoint, const unsigned int shape_nr);
231 
235  const Tensor<1, dim> &
236  derivative(const unsigned int qpoint, const unsigned int shape_nr) const;
237 
242  derivative(const unsigned int qpoint, const unsigned int shape_nr);
243 
247  const Tensor<2, dim> &
248  second_derivative(const unsigned int qpoint,
249  const unsigned int shape_nr) const;
250 
255  second_derivative(const unsigned int qpoint, const unsigned int shape_nr);
256 
260  const Tensor<3, dim> &
261  third_derivative(const unsigned int qpoint,
262  const unsigned int shape_nr) const;
263 
268  third_derivative(const unsigned int qpoint, const unsigned int shape_nr);
269 
273  const Tensor<4, dim> &
274  fourth_derivative(const unsigned int qpoint,
275  const unsigned int shape_nr) const;
276 
281  fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr);
282 
286  virtual std::size_t
287  memory_consumption() const override;
288 
294  std::vector<double> shape_values;
295 
301  std::vector<Tensor<1, dim>> shape_derivatives;
302 
309  std::vector<Tensor<2, dim>> shape_second_derivatives;
310 
317  std::vector<Tensor<3, dim>> shape_third_derivatives;
318 
325  std::vector<Tensor<4, dim>> shape_fourth_derivatives;
326 
333  std::array<std::vector<Tensor<1, dim>>,
336 
341 
345  const unsigned int polynomial_degree;
346 
350  const unsigned int n_shape_functions;
351 
361  mutable std::vector<DerivativeForm<1, dim, spacedim>> covariant;
362 
370  mutable std::vector<DerivativeForm<1, dim, spacedim>> contravariant;
371 
375  mutable std::vector<std::vector<Tensor<1, spacedim>>> aux;
376 
381  mutable std::vector<Point<spacedim>> mapping_support_points;
382 
388 
393  mutable std::vector<double> volume_elements;
394 
398  mutable std::vector<double> quadrature_weights;
399  };
400 
401 
402  // documentation can be found in Mapping::requires_update_flags()
403  virtual UpdateFlags
404  requires_update_flags(const UpdateFlags update_flags) const override;
405 
406  // documentation can be found in Mapping::get_data()
407  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
408  get_data(const UpdateFlags, const Quadrature<dim> &quadrature) const override;
409 
411 
412  // documentation can be found in Mapping::get_face_data()
413  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
414  get_face_data(const UpdateFlags flags,
415  const hp::QCollection<dim - 1> &quadrature) const override;
416 
417  // documentation can be found in Mapping::get_subface_data()
418  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
419  get_subface_data(const UpdateFlags flags,
420  const Quadrature<dim - 1> &quadrature) const override;
421 
422  // documentation can be found in Mapping::fill_fe_values()
425  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
426  const CellSimilarity::Similarity cell_similarity,
427  const Quadrature<dim> & quadrature,
428  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
430  &output_data) const override;
431 
433 
434  // documentation can be found in Mapping::fill_fe_face_values()
435  virtual void
437  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
438  const unsigned int face_no,
439  const hp::QCollection<dim - 1> & quadrature,
440  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
442  &output_data) const override;
443 
444  // documentation can be found in Mapping::fill_fe_subface_values()
445  virtual void
447  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
448  const unsigned int face_no,
449  const unsigned int subface_no,
450  const Quadrature<dim - 1> & quadrature,
451  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
453  &output_data) const override;
454 
459 protected:
460  const std::unique_ptr<FiniteElement<dim, spacedim>> fe;
461 
466  const unsigned int polynomial_degree;
467 
471  virtual std::vector<Point<spacedim>>
473  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const;
474 
475 private:
477 };
478 
479 
480 
483 /*----------------------------------------------------------------------*/
484 
485 #ifndef DOXYGEN
486 
487 template <int dim, int spacedim>
488 inline const double &
489 MappingFE<dim, spacedim>::InternalData::shape(const unsigned int qpoint,
490  const unsigned int shape_nr) const
491 {
492  AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
493  return shape_values[qpoint * n_shape_functions + shape_nr];
494 }
495 
496 
497 
498 template <int dim, int spacedim>
499 inline double &
500 MappingFE<dim, spacedim>::InternalData::shape(const unsigned int qpoint,
501  const unsigned int shape_nr)
502 {
503  AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
504  return shape_values[qpoint * n_shape_functions + shape_nr];
505 }
506 
507 
508 template <int dim, int spacedim>
509 inline const Tensor<1, dim> &
511  const unsigned int qpoint,
512  const unsigned int shape_nr) const
513 {
514  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
515  shape_derivatives.size());
516  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
517 }
518 
519 
520 
521 template <int dim, int spacedim>
522 inline Tensor<1, dim> &
524  const unsigned int shape_nr)
525 {
526  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
527  shape_derivatives.size());
528  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
529 }
530 
531 
532 template <int dim, int spacedim>
533 inline const Tensor<2, dim> &
535  const unsigned int qpoint,
536  const unsigned int shape_nr) const
537 {
538  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
539  shape_second_derivatives.size());
540  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
541 }
542 
543 
544 template <int dim, int spacedim>
545 inline Tensor<2, dim> &
547  const unsigned int qpoint,
548  const unsigned int shape_nr)
549 {
550  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
551  shape_second_derivatives.size());
552  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
553 }
554 
555 template <int dim, int spacedim>
556 inline const Tensor<3, dim> &
558  const unsigned int qpoint,
559  const unsigned int shape_nr) const
560 {
561  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
562  shape_third_derivatives.size());
563  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
564 }
565 
566 
567 template <int dim, int spacedim>
568 inline Tensor<3, dim> &
570  const unsigned int qpoint,
571  const unsigned int shape_nr)
572 {
573  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
574  shape_third_derivatives.size());
575  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
576 }
577 
578 
579 template <int dim, int spacedim>
580 inline const Tensor<4, dim> &
582  const unsigned int qpoint,
583  const unsigned int shape_nr) const
584 {
585  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
586  shape_fourth_derivatives.size());
587  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
588 }
589 
590 
591 template <int dim, int spacedim>
592 inline Tensor<4, dim> &
594  const unsigned int qpoint,
595  const unsigned int shape_nr)
596 {
597  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
598  shape_fourth_derivatives.size());
599  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
600 }
601 
602 
603 
604 template <int dim, int spacedim>
605 inline bool
607 {
608  return true;
609 }
610 
611 
612 
613 #endif // DOXYGEN
614 
615 /* -------------- declaration of explicit specializations ------------- */
616 
617 
619 
620 #endif
Triangulation< dim, spacedim >::cell_iterator cell_of_current_support_points
Definition: mapping_fe.h:387
void initialize_face(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
Definition: mapping_fe.cc:140
std::vector< Point< spacedim > > mapping_support_points
Definition: mapping_fe.h:381
const Tensor< 1, dim > & derivative(const unsigned int qpoint, const unsigned int shape_nr) const
Tensor< 2, dim > & second_derivative(const unsigned int qpoint, const unsigned int shape_nr)
std::vector< std::vector< Tensor< 1, spacedim > > > aux
Definition: mapping_fe.h:375
std::vector< Tensor< 3, dim > > shape_third_derivatives
Definition: mapping_fe.h:317
std::vector< Tensor< 2, dim > > shape_second_derivatives
Definition: mapping_fe.h:309
const Tensor< 3, dim > & third_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
std::vector< Tensor< 4, dim > > shape_fourth_derivatives
Definition: mapping_fe.h:325
void initialize(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
Definition: mapping_fe.cc:82
Tensor< 1, dim > & derivative(const unsigned int qpoint, const unsigned int shape_nr)
Tensor< 3, dim > & third_derivative(const unsigned int qpoint, const unsigned int shape_nr)
std::array< std::vector< Tensor< 1, dim > >, GeometryInfo< dim >::faces_per_cell *(dim - 1)> unit_tangentials
Definition: mapping_fe.h:335
std::vector< double > quadrature_weights
Definition: mapping_fe.h:398
const double & shape(const unsigned int qpoint, const unsigned int shape_nr) const
const FiniteElement< dim, spacedim > & fe
Definition: mapping_fe.h:340
std::vector< DerivativeForm< 1, dim, spacedim > > contravariant
Definition: mapping_fe.h:370
std::vector< DerivativeForm< 1, dim, spacedim > > covariant
Definition: mapping_fe.h:361
const Tensor< 4, dim > & fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
void compute_shape_function_values(const std::vector< Point< dim >> &unit_points)
Definition: mapping_fe.cc:180
std::vector< double > shape_values
Definition: mapping_fe.h:294
std::vector< Tensor< 1, dim > > shape_derivatives
Definition: mapping_fe.h:301
Tensor< 4, dim > & fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr)
virtual std::size_t memory_consumption() const override
Definition: mapping_fe.cc:62
double & shape(const unsigned int qpoint, const unsigned int shape_nr)
InternalData(const FiniteElement< dim, spacedim > &fe)
Definition: mapping_fe.cc:51
const unsigned int n_shape_functions
Definition: mapping_fe.h:350
std::vector< double > volume_elements
Definition: mapping_fe.h:393
const Tensor< 2, dim > & second_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
const unsigned int polynomial_degree
Definition: mapping_fe.h:345
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_fe.cc:1635
const unsigned int polynomial_degree
Definition: mapping_fe.h:466
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_fe.cc:1589
MappingFE(const FiniteElement< dim, spacedim > &fe)
Definition: mapping_fe.cc:845
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_fe.cc:1075
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const override
Definition: mapping_fe.cc:923
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
Definition: mapping_fe.cc:1004
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_fe.cc:2091
Table< 2, double > mapping_support_point_weights
Definition: mapping_fe.h:476
virtual std::vector< Point< spacedim > > compute_mapping_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const
Definition: mapping_fe.cc:2278
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_fe.cc:1113
virtual BoundingBox< spacedim > get_bounding_box(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const override
Definition: mapping_fe.cc:2305
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const override
Definition: mapping_fe.cc:1060
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags flags, const Quadrature< dim - 1 > &quadrature) const override
Definition: mapping_fe.cc:1094
unsigned int get_degree() const
Definition: mapping_fe.cc:895
virtual bool is_compatible_with(const ReferenceCell &reference_cell) const override
Definition: mapping_fe.cc:2315
virtual bool preserves_vertex_locations() const override
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
Definition: mapping_fe.cc:886
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const override
Definition: mapping_fe.cc:904
const std::unique_ptr< FiniteElement< dim, spacedim > > fe
Definition: mapping_fe.h:460
Abstract base class for mapping classes.
Definition: mapping.h:311
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:458
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:459
#define AssertIndexRange(index, range)
Definition: exceptions.h:1760
UpdateFlags
MappingKind
Definition: mapping.h:72
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)