Reference documentation for deal.II version Git 08727cc441 2020-07-02 15:45:42 -0400
\(\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_field.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2001 - 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 
26 
27 #include <deal.II/fe/fe.h>
28 #include <deal.II/fe/fe_values.h>
29 #include <deal.II/fe/mapping.h>
30 
31 #include <deal.II/lac/vector.h>
32 
33 #include <array>
34 
35 
37 
38 
41 
77 template <int dim,
78  int spacedim = dim,
79  typename VectorType = Vector<double>,
80  typename DoFHandlerType = DoFHandler<dim, spacedim>>
81 class MappingFEField : public Mapping<dim, spacedim>
82 {
83 public:
116  MappingFEField(const DoFHandlerType &euler_dof_handler,
117  const VectorType & euler_vector,
118  const ComponentMask & mask = ComponentMask());
119 
129  MappingFEField(const DoFHandlerType & euler_dof_handler,
130  const std::vector<VectorType> &euler_vector,
131  const ComponentMask & mask = ComponentMask());
132 
140  MappingFEField(const DoFHandlerType & euler_dof_handler,
141  const MGLevelObject<VectorType> &euler_vector,
142  const ComponentMask & mask = ComponentMask());
143 
149 
154  virtual std::unique_ptr<Mapping<dim, spacedim>>
155  clone() const override;
156 
162  virtual bool
163  preserves_vertex_locations() const override;
164 
172  virtual std::array<Point<spacedim>, GeometryInfo<dim>::vertices_per_cell>
174  const override;
175 
181  // for documentation, see the Mapping base class
182  virtual Point<spacedim>
184  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
185  const Point<dim> &p) const override;
186 
187  // for documentation, see the Mapping base class
188  virtual Point<dim>
190  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
191  const Point<spacedim> &p) const override;
192 
202  // for documentation, see the Mapping base class
203  virtual void
204  transform(const ArrayView<const Tensor<1, dim>> & input,
205  const MappingKind kind,
206  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
207  const ArrayView<Tensor<1, spacedim>> &output) const override;
208 
209  // for documentation, see the Mapping base class
210  virtual void
212  const MappingKind kind,
213  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
214  const ArrayView<Tensor<2, spacedim>> &output) const override;
215 
216  // for documentation, see the Mapping base class
217  virtual void
218  transform(const ArrayView<const Tensor<2, dim>> & input,
219  const MappingKind kind,
220  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
221  const ArrayView<Tensor<2, spacedim>> &output) const override;
222 
223  // for documentation, see the Mapping base class
224  virtual void
226  const MappingKind kind,
227  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
228  const ArrayView<Tensor<3, spacedim>> &output) const override;
229 
230  // for documentation, see the Mapping base class
231  virtual void
232  transform(const ArrayView<const Tensor<3, dim>> & input,
233  const MappingKind kind,
234  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
235  const ArrayView<Tensor<3, spacedim>> &output) const override;
236 
245  unsigned int
246  get_degree() const;
247 
253  get_component_mask() const;
254 
259 
260 private:
266  // documentation can be found in Mapping::requires_update_flags()
267  virtual UpdateFlags
268  requires_update_flags(const UpdateFlags update_flags) const override;
269 
270 public:
282  class InternalData : public Mapping<dim, spacedim>::InternalDataBase
283  {
284  public:
289  const ComponentMask & mask);
290 
295  const double &
296  shape(const unsigned int qpoint, const unsigned int shape_nr) const;
297 
301  double &
302  shape(const unsigned int qpoint, const unsigned int shape_nr);
303 
307  const Tensor<1, dim> &
308  derivative(const unsigned int qpoint, const unsigned int shape_nr) const;
309 
314  derivative(const unsigned int qpoint, const unsigned int shape_nr);
315 
319  const Tensor<2, dim> &
320  second_derivative(const unsigned int qpoint,
321  const unsigned int shape_nr) const;
322 
327  second_derivative(const unsigned int qpoint, const unsigned int shape_nr);
328 
332  const Tensor<3, dim> &
333  third_derivative(const unsigned int qpoint,
334  const unsigned int shape_nr) const;
335 
340  third_derivative(const unsigned int qpoint, const unsigned int shape_nr);
341 
345  const Tensor<4, dim> &
346  fourth_derivative(const unsigned int qpoint,
347  const unsigned int shape_nr) const;
348 
353  fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr);
354 
358  virtual std::size_t
359  memory_consumption() const override;
360 
366  std::vector<double> shape_values;
367 
373  std::vector<Tensor<1, dim>> shape_derivatives;
374 
381  std::vector<Tensor<2, dim>> shape_second_derivatives;
382 
389  std::vector<Tensor<3, dim>> shape_third_derivatives;
390 
397  std::vector<Tensor<4, dim>> shape_fourth_derivatives;
398 
412  std::array<std::vector<Tensor<1, dim>>,
415 
422  unsigned int n_shape_functions;
423 
437 
447  mutable std::vector<DerivativeForm<1, dim, spacedim>> covariant;
448 
456  mutable std::vector<DerivativeForm<1, dim, spacedim>> contravariant;
457 
462  mutable std::vector<double> volume_elements;
463 
467  mutable std::vector<std::vector<Tensor<1, spacedim>>> aux;
468 
472  mutable std::vector<types::global_dof_index> local_dof_indices;
473 
477  mutable std::vector<double> local_dof_values;
478  };
479 
480 private:
481  // documentation can be found in Mapping::get_data()
482  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
483  get_data(const UpdateFlags, const Quadrature<dim> &quadrature) const override;
484 
485  // documentation can be found in Mapping::get_face_data()
486  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
487  get_face_data(const UpdateFlags flags,
488  const Quadrature<dim - 1> &quadrature) const override;
489 
490  // documentation can be found in Mapping::get_subface_data()
491  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
492  get_subface_data(const UpdateFlags flags,
493  const Quadrature<dim - 1> &quadrature) const override;
494 
495  // documentation can be found in Mapping::fill_fe_values()
498  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
499  const CellSimilarity::Similarity cell_similarity,
500  const Quadrature<dim> & quadrature,
501  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
503  &output_data) const override;
504 
505  // documentation can be found in Mapping::fill_fe_face_values()
506  virtual void
508  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
509  const unsigned int face_no,
510  const Quadrature<dim - 1> & quadrature,
511  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
513  &output_data) const override;
514 
515  // documentation can be found in Mapping::fill_fe_subface_values()
516  virtual void
518  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
519  const unsigned int face_no,
520  const unsigned int subface_no,
521  const Quadrature<dim - 1> & quadrature,
522  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
524  &output_data) const override;
525 
534  const bool uses_level_dofs;
535 
539  std::vector<
540  SmartPointer<const VectorType,
543 
547  SmartPointer<const DoFHandlerType,
550 
551 private:
568  do_transform_unit_to_real_cell(const InternalData &mdata) const;
569 
570 
585  Point<dim>
587  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
588  const Point<spacedim> & p,
589  const Point<dim> & initial_p_unit,
590  InternalData & mdata) const;
591 
595  void
597  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
599  InternalData &data) const;
600 
604  virtual void
606  const std::vector<Point<dim>> &unit_points,
608  InternalData &data) const;
609 
610  /*
611  * Which components to use for the mapping.
612  */
614 
624  std::vector<unsigned int> fe_to_real;
625 
631 
635  mutable std::mutex fe_values_mutex;
636 
637  void
638  compute_data(const UpdateFlags update_flags,
639  const Quadrature<dim> &q,
640  const unsigned int n_original_q_points,
641  InternalData & data) const;
642 
643  void
644  compute_face_data(const UpdateFlags update_flags,
645  const Quadrature<dim> &q,
646  const unsigned int n_original_q_points,
647  InternalData & data) const;
648 
649 
650  // Declare other MappingFEField classes friends.
651  template <int, int, class, class>
652  friend class MappingFEField;
653 };
654 
657 /* -------------- declaration of explicit specializations ------------- */
658 
659 #ifndef DOXYGEN
660 
661 #endif // DOXYGEN
662 
664 
665 #endif
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const override
void compute_face_data(const UpdateFlags update_flags, const Quadrature< dim > &q, const unsigned int n_original_q_points, InternalData &data) const
std::mutex fe_values_mutex
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
virtual void compute_shapes_virtual(const std::vector< Point< dim >> &unit_points, typename MappingFEField< dim, spacedim, VectorType, DoFHandlerType >::InternalData &data) const
std::vector< unsigned int > fe_to_real
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
const ComponentMask fe_mask
MappingKind
Definition: mapping.h:62
Point< spacedim > do_transform_unit_to_real_cell(const InternalData &mdata) const
unsigned int get_degree() const
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< 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
std::vector< Tensor< 3, dim > > shape_third_derivatives
std::vector< Tensor< 1, dim > > shape_derivatives
const double & shape(const unsigned int qpoint, const unsigned int shape_nr) const
std::vector< Tensor< 2, dim > > shape_second_derivatives
UpdateFlags
std::vector< Tensor< 4, dim > > shape_fourth_derivatives
Abstract base class for mapping classes.
Definition: mapping.h:301
#define DeclException0(Exception0)
Definition: exceptions.h:470
std::vector< double > volume_elements
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_face_data(const UpdateFlags flags, const Quadrature< dim - 1 > &quadrature) const override
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const override
std::vector< DerivativeForm< 1, dim, spacedim > > covariant
std::vector< double > local_dof_values
const Tensor< 3, dim > & third_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
virtual std::array< Point< spacedim >, GeometryInfo< dim >::vertices_per_cell > get_vertices(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const override
std::vector< DerivativeForm< 1, dim, spacedim > > contravariant
std::vector< double > shape_values
std::vector< std::vector< Tensor< 1, spacedim > > > aux
friend class MappingFEField
std::vector< types::global_dof_index > local_dof_indices
void compute_data(const UpdateFlags update_flags, const Quadrature< dim > &q, const unsigned int n_original_q_points, InternalData &data) const
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
const Tensor< 2, dim > & second_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
const bool uses_level_dofs
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:362
const Tensor< 4, dim > & fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
SmartPointer< const DoFHandlerType, MappingFEField< dim, spacedim, VectorType, DoFHandlerType > > euler_dof_handler
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, DoFHandlerType > > > euler_vector
const Tensor< 1, dim > & derivative(const unsigned int qpoint, const unsigned int shape_nr) const
ComponentMask get_component_mask() const
void update_internal_dofs(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const typename MappingFEField< dim, spacedim, VectorType, DoFHandlerType >::InternalData &data) const
virtual bool preserves_vertex_locations() const override
static ::ExceptionBase & ExcInactiveCell()
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
std::array< std::vector< Tensor< 1, dim > >, GeometryInfo< dim >::faces_per_cell *(dim - 1)> unit_tangentials
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 std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags flags, const Quadrature< dim - 1 > &quadrature) const override
virtual std::size_t memory_consumption() const override
InternalData(const FiniteElement< dim, spacedim > &fe, const ComponentMask &mask)
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override