Reference documentation for deal.II version Git 354f154960 2021-01-18 17:44:59 +0100
\(\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 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 
39 
40 
57 template <int dim, int spacedim = dim>
58 class MappingFE : public Mapping<dim, spacedim>
59 {
60 public:
64  explicit MappingFE(const FiniteElement<dim, spacedim> &fe);
65 
69  MappingFE(const MappingFE<dim, spacedim> &mapping);
70 
71  // for documentation, see the Mapping base class
72  virtual std::unique_ptr<Mapping<dim, spacedim>>
73  clone() const override;
74 
79  unsigned int
80  get_degree() const;
81 
82  // for documentation, see the Mapping base class
83  virtual BoundingBox<spacedim>
85  &cell) const override;
86 
87 
92  virtual bool
93  preserves_vertex_locations() const override;
94 
100  // for documentation, see the Mapping base class
101  virtual Point<spacedim>
103  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
104  const Point<dim> &p) const override;
105 
106  // for documentation, see the Mapping base class
107  virtual Point<dim>
109  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
110  const Point<spacedim> &p) const override;
111 
121  // for documentation, see the Mapping base class
122  virtual void
123  transform(const ArrayView<const Tensor<1, dim>> & input,
124  const MappingKind kind,
125  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
126  const ArrayView<Tensor<1, spacedim>> &output) const override;
127 
128  // for documentation, see the Mapping base class
129  virtual void
131  const MappingKind kind,
132  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
133  const ArrayView<Tensor<2, spacedim>> &output) const override;
134 
135  // for documentation, see the Mapping base class
136  virtual void
137  transform(const ArrayView<const Tensor<2, dim>> & input,
138  const MappingKind kind,
139  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
140  const ArrayView<Tensor<2, spacedim>> &output) const override;
141 
142  // for documentation, see the Mapping base class
143  virtual void
145  const MappingKind kind,
146  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
147  const ArrayView<Tensor<3, spacedim>> &output) const override;
148 
149  // for documentation, see the Mapping base class
150  virtual void
151  transform(const ArrayView<const Tensor<3, dim>> & input,
152  const MappingKind kind,
153  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
154  const ArrayView<Tensor<3, spacedim>> &output) const override;
155 
176  class InternalData : public Mapping<dim, spacedim>::InternalDataBase
177  {
178  public:
183 
192  void
193  initialize(const UpdateFlags update_flags,
194  const Quadrature<dim> &quadrature,
195  const unsigned int n_original_q_points);
196 
202  void
203  initialize_face(const UpdateFlags update_flags,
204  const Quadrature<dim> &quadrature,
205  const unsigned int n_original_q_points);
206 
211  void
212  compute_shape_function_values(const std::vector<Point<dim>> &unit_points);
213 
214 
219  const double &
220  shape(const unsigned int qpoint, const unsigned int shape_nr) const;
221 
225  double &
226  shape(const unsigned int qpoint, const unsigned int shape_nr);
227 
231  const Tensor<1, dim> &
232  derivative(const unsigned int qpoint, const unsigned int shape_nr) const;
233 
238  derivative(const unsigned int qpoint, const unsigned int shape_nr);
239 
243  const Tensor<2, dim> &
244  second_derivative(const unsigned int qpoint,
245  const unsigned int shape_nr) const;
246 
251  second_derivative(const unsigned int qpoint, const unsigned int shape_nr);
252 
256  const Tensor<3, dim> &
257  third_derivative(const unsigned int qpoint,
258  const unsigned int shape_nr) const;
259 
264  third_derivative(const unsigned int qpoint, const unsigned int shape_nr);
265 
269  const Tensor<4, dim> &
270  fourth_derivative(const unsigned int qpoint,
271  const unsigned int shape_nr) const;
272 
277  fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr);
278 
282  virtual std::size_t
283  memory_consumption() const override;
284 
290  std::vector<double> shape_values;
291 
297  std::vector<Tensor<1, dim>> shape_derivatives;
298 
305  std::vector<Tensor<2, dim>> shape_second_derivatives;
306 
313  std::vector<Tensor<3, dim>> shape_third_derivatives;
314 
321  std::vector<Tensor<4, dim>> shape_fourth_derivatives;
322 
329  std::array<std::vector<Tensor<1, dim>>,
332 
337 
341  const unsigned int polynomial_degree;
342 
346  const unsigned int n_shape_functions;
347 
357  mutable std::vector<DerivativeForm<1, dim, spacedim>> covariant;
358 
366  mutable std::vector<DerivativeForm<1, dim, spacedim>> contravariant;
367 
371  mutable std::vector<std::vector<Tensor<1, spacedim>>> aux;
372 
377  mutable std::vector<Point<spacedim>> mapping_support_points;
378 
384 
389  mutable std::vector<double> volume_elements;
390 
394  mutable std::vector<double> quadrature_weights;
395  };
396 
397 
398  // documentation can be found in Mapping::requires_update_flags()
399  virtual UpdateFlags
400  requires_update_flags(const UpdateFlags update_flags) const override;
401 
402  // documentation can be found in Mapping::get_data()
403  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
404  get_data(const UpdateFlags, const Quadrature<dim> &quadrature) const override;
405 
407 
408  // documentation can be found in Mapping::get_face_data()
409  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
410  get_face_data(const UpdateFlags flags,
411  const hp::QCollection<dim - 1> &quadrature) const override;
412 
413  // documentation can be found in Mapping::get_subface_data()
414  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
415  get_subface_data(const UpdateFlags flags,
416  const Quadrature<dim - 1> &quadrature) const override;
417 
418  // documentation can be found in Mapping::fill_fe_values()
421  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
422  const CellSimilarity::Similarity cell_similarity,
423  const Quadrature<dim> & quadrature,
424  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
426  &output_data) const override;
427 
429 
430  // documentation can be found in Mapping::fill_fe_face_values()
431  virtual void
433  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
434  const unsigned int face_no,
435  const hp::QCollection<dim - 1> & quadrature,
436  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
438  &output_data) const override;
439 
440  // documentation can be found in Mapping::fill_fe_subface_values()
441  virtual void
443  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
444  const unsigned int face_no,
445  const unsigned int subface_no,
446  const Quadrature<dim - 1> & quadrature,
447  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
449  &output_data) const override;
450 
455 protected:
456  const std::unique_ptr<FiniteElement<dim, spacedim>> fe;
457 
462  const unsigned int polynomial_degree;
463 
467  virtual std::vector<Point<spacedim>>
469  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const;
470 
471 private:
473 };
474 
475 
476 
479 /*----------------------------------------------------------------------*/
480 
481 #ifndef DOXYGEN
482 
483 template <int dim, int spacedim>
484 inline const double &
485 MappingFE<dim, spacedim>::InternalData::shape(const unsigned int qpoint,
486  const unsigned int shape_nr) const
487 {
488  AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
489  return shape_values[qpoint * n_shape_functions + shape_nr];
490 }
491 
492 
493 
494 template <int dim, int spacedim>
495 inline double &
496 MappingFE<dim, spacedim>::InternalData::shape(const unsigned int qpoint,
497  const unsigned int shape_nr)
498 {
499  AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
500  return shape_values[qpoint * n_shape_functions + shape_nr];
501 }
502 
503 
504 template <int dim, int spacedim>
505 inline const Tensor<1, dim> &
507  const unsigned int qpoint,
508  const unsigned int shape_nr) const
509 {
510  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
511  shape_derivatives.size());
512  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
513 }
514 
515 
516 
517 template <int dim, int spacedim>
518 inline Tensor<1, dim> &
520  const unsigned int shape_nr)
521 {
522  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
523  shape_derivatives.size());
524  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
525 }
526 
527 
528 template <int dim, int spacedim>
529 inline const Tensor<2, dim> &
531  const unsigned int qpoint,
532  const unsigned int shape_nr) const
533 {
534  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
535  shape_second_derivatives.size());
536  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
537 }
538 
539 
540 template <int dim, int spacedim>
541 inline Tensor<2, dim> &
543  const unsigned int qpoint,
544  const unsigned int shape_nr)
545 {
546  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
547  shape_second_derivatives.size());
548  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
549 }
550 
551 template <int dim, int spacedim>
552 inline const Tensor<3, dim> &
554  const unsigned int qpoint,
555  const unsigned int shape_nr) const
556 {
557  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
558  shape_third_derivatives.size());
559  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
560 }
561 
562 
563 template <int dim, int spacedim>
564 inline Tensor<3, dim> &
566  const unsigned int qpoint,
567  const unsigned int shape_nr)
568 {
569  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
570  shape_third_derivatives.size());
571  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
572 }
573 
574 
575 template <int dim, int spacedim>
576 inline const Tensor<4, dim> &
578  const unsigned int qpoint,
579  const unsigned int shape_nr) const
580 {
581  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
582  shape_fourth_derivatives.size());
583  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
584 }
585 
586 
587 template <int dim, int spacedim>
588 inline Tensor<4, dim> &
590  const unsigned int qpoint,
591  const unsigned int shape_nr)
592 {
593  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
594  shape_fourth_derivatives.size());
595  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
596 }
597 
598 
599 
600 template <int dim, int spacedim>
601 inline bool
603 {
604  return true;
605 }
606 
607 #endif // DOXYGEN
608 
609 /* -------------- declaration of explicit specializations ------------- */
610 
611 
613 
614 #endif
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
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const override
Definition: mapping_fe.cc:1066
std::vector< Tensor< 4, dim > > shape_fourth_derivatives
Definition: mapping_fe.h:321
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
Definition: mapping_fe.cc:892
MappingFE(const FiniteElement< dim, spacedim > &fe)
Definition: mapping_fe.cc:849
std::vector< double > volume_elements
Definition: mapping_fe.h:389
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:1119
std::vector< std::vector< Tensor< 1, spacedim > > > aux
Definition: mapping_fe.h:371
std::vector< double > shape_values
Definition: mapping_fe.h:290
InternalData(const FiniteElement< dim, spacedim > &fe)
Definition: mapping_fe.cc:50
#define AssertIndexRange(index, range)
Definition: exceptions.h:1691
unsigned int get_degree() const
Definition: mapping_fe.cc:901
virtual BoundingBox< spacedim > get_bounding_box(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const override
Definition: mapping_fe.cc:2305
virtual bool preserves_vertex_locations() const override
MappingKind
Definition: mapping.h:64
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:929
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:1641
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:910
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:1100
std::vector< Tensor< 3, dim > > shape_third_derivatives
Definition: mapping_fe.h:313
void initialize_face(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
Definition: mapping_fe.cc:139
std::vector< Tensor< 1, dim > > shape_derivatives
Definition: mapping_fe.h:297
UpdateFlags
void initialize(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
Definition: mapping_fe.cc:81
Abstract base class for mapping classes.
Definition: mapping.h:303
virtual std::size_t memory_consumption() const override
Definition: mapping_fe.cc:61
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:380
std::vector< Tensor< 2, dim > > shape_second_derivatives
Definition: mapping_fe.h:305
Triangulation< dim, spacedim >::cell_iterator cell_of_current_support_points
Definition: mapping_fe.h:383
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:1595
std::vector< Point< spacedim > > mapping_support_points
Definition: mapping_fe.h:377
const double & shape(const unsigned int qpoint, const unsigned int shape_nr) const
std::vector< DerivativeForm< 1, dim, spacedim > > covariant
Definition: mapping_fe.h:357
const std::unique_ptr< FiniteElement< dim, spacedim > > fe
Definition: mapping_fe.h:456
const Tensor< 3, dim > & third_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
virtual std::vector< Point< spacedim > > compute_mapping_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const
Definition: mapping_fe.cc:2278
const Tensor< 2, dim > & second_derivative(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
Definition: mapping_fe.h:331
const Tensor< 1, dim > & derivative(const unsigned int qpoint, const unsigned int shape_nr) const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:379
void compute_shape_function_values(const std::vector< Point< dim >> &unit_points)
Definition: mapping_fe.cc:181
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:1081
const unsigned int polynomial_degree
Definition: mapping_fe.h:341
std::vector< double > quadrature_weights
Definition: mapping_fe.h:394
Table< 2, double > mapping_support_point_weights
Definition: mapping_fe.h:472
const unsigned int n_shape_functions
Definition: mapping_fe.h:346
const unsigned int polynomial_degree
Definition: mapping_fe.h:462
const Tensor< 4, dim > & fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
std::vector< DerivativeForm< 1, dim, spacedim > > contravariant
Definition: mapping_fe.h:366
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
Definition: mapping_fe.cc:1010
const FiniteElement< dim, spacedim > & fe
Definition: mapping_fe.h:336