Reference documentation for deal.II version 9.5.0
\(\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\}}\)
Loading...
Searching...
No Matches
mapping_fe_field.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2001 - 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_fe_field_h
17#define dealii_mapping_fe_field_h
18
19
20#include <deal.II/base/config.h>
21
23#include <deal.II/base/mutex.h>
24
26
27#include <deal.II/fe/fe.h>
29#include <deal.II/fe/mapping.h>
30
31#include <deal.II/lac/vector.h>
32
33#include <array>
34
35
37
38
79template <int dim, int spacedim = dim, typename VectorType = Vector<double>>
80class MappingFEField : public Mapping<dim, spacedim>
81{
82public:
116 const VectorType & euler_vector,
117 const ComponentMask & mask = ComponentMask());
118
129 const std::vector<VectorType> & euler_vector,
130 const ComponentMask & mask = ComponentMask());
131
141 const ComponentMask & mask = ComponentMask());
142
147
152 virtual std::unique_ptr<Mapping<dim, spacedim>>
153 clone() const override;
154
160 virtual bool
161 preserves_vertex_locations() const override;
162
163 virtual bool
164 is_compatible_with(const ReferenceCell &reference_cell) const override;
165
173 virtual boost::container::small_vector<Point<spacedim>,
176 const override;
177
183 // for documentation, see the Mapping base class
184 virtual Point<spacedim>
187 const Point<dim> &p) const override;
188
189 // for documentation, see the Mapping base class
190 virtual Point<dim>
193 const Point<spacedim> &p) const override;
194
204 // for documentation, see the Mapping base class
205 virtual void
206 transform(const ArrayView<const Tensor<1, dim>> & input,
207 const MappingKind kind,
209 const ArrayView<Tensor<1, spacedim>> &output) const override;
210
211 // for documentation, see the Mapping base class
212 virtual void
214 const MappingKind kind,
216 const ArrayView<Tensor<2, spacedim>> &output) const override;
217
218 // for documentation, see the Mapping base class
219 virtual void
220 transform(const ArrayView<const Tensor<2, dim>> & input,
221 const MappingKind kind,
223 const ArrayView<Tensor<2, spacedim>> &output) const override;
224
225 // for documentation, see the Mapping base class
226 virtual void
228 const MappingKind kind,
230 const ArrayView<Tensor<3, spacedim>> &output) const override;
231
232 // for documentation, see the Mapping base class
233 virtual void
234 transform(const ArrayView<const Tensor<3, dim>> & input,
235 const MappingKind kind,
237 const ArrayView<Tensor<3, spacedim>> &output) const override;
238
247 unsigned int
248 get_degree() const;
249
255 get_component_mask() const;
256
261
262private:
268 // documentation can be found in Mapping::requires_update_flags()
269 virtual UpdateFlags
270 requires_update_flags(const UpdateFlags update_flags) const override;
271
272public:
284 class InternalData : public Mapping<dim, spacedim>::InternalDataBase
285 {
286 public:
291 const ComponentMask & mask);
292
297 const double &
298 shape(const unsigned int qpoint, const unsigned int shape_nr) const;
299
303 double &
304 shape(const unsigned int qpoint, const unsigned int shape_nr);
305
309 const Tensor<1, dim> &
310 derivative(const unsigned int qpoint, const unsigned int shape_nr) const;
311
316 derivative(const unsigned int qpoint, const unsigned int shape_nr);
317
321 const Tensor<2, dim> &
322 second_derivative(const unsigned int qpoint,
323 const unsigned int shape_nr) const;
324
329 second_derivative(const unsigned int qpoint, const unsigned int shape_nr);
330
334 const Tensor<3, dim> &
335 third_derivative(const unsigned int qpoint,
336 const unsigned int shape_nr) const;
337
342 third_derivative(const unsigned int qpoint, const unsigned int shape_nr);
343
347 const Tensor<4, dim> &
348 fourth_derivative(const unsigned int qpoint,
349 const unsigned int shape_nr) const;
350
355 fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr);
356
360 virtual std::size_t
361 memory_consumption() const override;
362
368 std::vector<double> shape_values;
369
375 std::vector<Tensor<1, dim>> shape_derivatives;
376
383 std::vector<Tensor<2, dim>> shape_second_derivatives;
384
391 std::vector<Tensor<3, dim>> shape_third_derivatives;
392
399 std::vector<Tensor<4, dim>> shape_fourth_derivatives;
400
414 std::array<std::vector<Tensor<1, dim>>,
417
424 unsigned int n_shape_functions;
425
439
449 mutable std::vector<DerivativeForm<1, dim, spacedim>> covariant;
450
458 mutable std::vector<DerivativeForm<1, dim, spacedim>> contravariant;
459
464 mutable std::vector<double> volume_elements;
465
469 mutable std::vector<double> quadrature_weights;
470
474 mutable std::vector<std::vector<Tensor<1, spacedim>>> aux;
475
479 mutable std::vector<types::global_dof_index> local_dof_indices;
480
484 mutable std::vector<double> local_dof_values;
485 };
486
487protected:
488 // documentation can be found in Mapping::get_data()
489 virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
490 get_data(const UpdateFlags, const Quadrature<dim> &quadrature) const override;
491
492 using Mapping<dim, spacedim>::get_face_data;
493
494 // documentation can be found in Mapping::get_face_data()
495 virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
496 get_face_data(const UpdateFlags flags,
497 const hp::QCollection<dim - 1> &quadrature) const override;
498
499 // documentation can be found in Mapping::get_subface_data()
500 virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
501 get_subface_data(const UpdateFlags flags,
502 const Quadrature<dim - 1> &quadrature) const override;
503
504 // documentation can be found in Mapping::fill_fe_values()
508 const CellSimilarity::Similarity cell_similarity,
509 const Quadrature<dim> & quadrature,
510 const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
512 &output_data) const override;
513
514 using Mapping<dim, spacedim>::fill_fe_face_values;
515
516 // documentation can be found in Mapping::fill_fe_face_values()
517 virtual void
520 const unsigned int face_no,
521 const hp::QCollection<dim - 1> & quadrature,
522 const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
524 &output_data) const override;
525
526 // documentation can be found in Mapping::fill_fe_subface_values()
527 virtual void
530 const unsigned int face_no,
531 const unsigned int subface_no,
532 const Quadrature<dim - 1> & quadrature,
533 const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
535 &output_data) const override;
536
537 virtual void
541 const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
543 &output_data) const override;
544
554
559 const bool uses_level_dofs;
560
564 std::vector<
567
574
575private:
593
594
612 const Point<spacedim> & p,
613 const Point<dim> & initial_p_unit,
614 InternalData & mdata) const;
615
619 void
623 &data) const;
624
628 virtual void
630 const std::vector<Point<dim>> &unit_points,
632 const;
633
634 /*
635 * Which components to use for the mapping.
636 */
638
648 std::vector<unsigned int> fe_to_real;
649
655
660
661 void
662 compute_data(const UpdateFlags update_flags,
663 const Quadrature<dim> &q,
664 const unsigned int n_original_q_points,
665 InternalData & data) const;
666
667 void
668 compute_face_data(const UpdateFlags update_flags,
669 const Quadrature<dim> &q,
670 const unsigned int n_original_q_points,
671 InternalData & data) const;
672
673
674 // Declare other MappingFEField classes friends.
675 template <int, int, class>
676 friend class MappingFEField;
677};
678
681/* -------------- declaration of explicit specializations ------------- */
682
683#ifndef DOXYGEN
684
685#endif // DOXYGEN
686
688
689#endif
const Tensor< 3, dim > & third_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
std::vector< double > quadrature_weights
std::vector< DerivativeForm< 1, dim, spacedim > > contravariant
std::vector< std::vector< Tensor< 1, spacedim > > > aux
std::vector< Tensor< 4, dim > > shape_fourth_derivatives
std::vector< Tensor< 3, dim > > shape_third_derivatives
const Tensor< 2, dim > & second_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
std::vector< types::global_dof_index > local_dof_indices
virtual std::size_t memory_consumption() const override
const Tensor< 4, dim > & fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
std::vector< Tensor< 1, dim > > shape_derivatives
std::vector< DerivativeForm< 1, dim, spacedim > > covariant
const Tensor< 1, dim > & derivative(const unsigned int qpoint, const unsigned int shape_nr) const
std::vector< double > local_dof_values
const double & shape(const unsigned int qpoint, const unsigned int shape_nr) const
std::vector< double > shape_values
std::vector< Tensor< 2, dim > > shape_second_derivatives
std::vector< double > volume_elements
virtual boost::container::small_vector< Point< spacedim >, GeometryInfo< dim >::vertices_per_cell > get_vertices(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const override
std::vector< unsigned int > fe_to_real
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags flags, const Quadrature< dim - 1 > &quadrature) const override
const ComponentMask fe_mask
void update_internal_dofs(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const typename MappingFEField< dim, spacedim, VectorType >::InternalData &data) const
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
void compute_face_data(const UpdateFlags update_flags, const Quadrature< dim > &q, const unsigned int n_original_q_points, InternalData &data) const
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const override
std::vector< SmartPointer< const VectorType, MappingFEField< dim, spacedim, VectorType > > > euler_vector
virtual bool preserves_vertex_locations() const override
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const override
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_face_data(const UpdateFlags flags, const hp::QCollection< dim - 1 > &quadrature) const override
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
friend class MappingFEField
virtual void compute_shapes_virtual(const std::vector< Point< dim > > &unit_points, typename MappingFEField< dim, spacedim, VectorType >::InternalData &data) const
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const override
virtual bool is_compatible_with(const ReferenceCell &reference_cell) const override
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
void compute_data(const UpdateFlags update_flags, const Quadrature< dim > &q, const unsigned int n_original_q_points, InternalData &data) const
ReferenceCell reference_cell
SmartPointer< const DoFHandler< dim, spacedim >, MappingFEField< dim, spacedim, VectorType > > euler_dof_handler
unsigned int get_degree() const
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
Point< spacedim > do_transform_unit_to_real_cell(const InternalData &mdata) const
Threads::Mutex fe_values_mutex
const bool uses_level_dofs
ComponentMask get_component_mask() const
Point< dim > do_transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p, const Point< dim > &initial_p_unit, InternalData &mdata) const
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
FEValues< dim, spacedim > fe_values
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
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
Abstract base class for mapping classes.
Definition mapping.h:317
Definition point.h:112
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
#define DeclException0(Exception0)
Definition exceptions.h:465
static ::ExceptionBase & ExcInactiveCell()
UpdateFlags
MappingKind
Definition mapping.h:78