Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
mapping_fe_field.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2001 - 2018 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 
22 #include <deal.II/base/thread_management.h>
23 
24 #include <deal.II/dofs/dof_handler.h>
25 
26 #include <deal.II/fe/fe.h>
27 #include <deal.II/fe/fe_values.h>
28 #include <deal.II/fe/mapping.h>
29 
30 #include <deal.II/lac/vector.h>
31 
32 #include <array>
33 
34 
35 DEAL_II_NAMESPACE_OPEN
36 
37 
40 
78 template <int dim,
79  int spacedim = dim,
80  typename VectorType = Vector<double>,
81  typename DoFHandlerType = DoFHandler<dim, spacedim>>
82 class MappingFEField : public Mapping<dim, spacedim>
83 {
84 public:
117  MappingFEField(const DoFHandlerType &euler_dof_handler,
118  const VectorType & euler_vector,
119  const ComponentMask & mask = ComponentMask());
120 
126 
131  virtual std::unique_ptr<Mapping<dim, spacedim>>
132  clone() const override;
133 
139  virtual bool
140  preserves_vertex_locations() const override;
141 
149  virtual std::array<Point<spacedim>, GeometryInfo<dim>::vertices_per_cell>
151  const override;
152 
158  // for documentation, see the Mapping base class
159  virtual Point<spacedim>
161  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
162  const Point<dim> &p) const override;
163 
164  // for documentation, see the Mapping base class
165  virtual Point<dim>
167  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
168  const Point<spacedim> &p) const override;
169 
179  // for documentation, see the Mapping base class
180  virtual void
181  transform(const ArrayView<const Tensor<1, dim>> & input,
182  const MappingType type,
183  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
184  const ArrayView<Tensor<1, spacedim>> &output) const override;
185 
186  // for documentation, see the Mapping base class
187  virtual void
189  const MappingType type,
190  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
191  const ArrayView<Tensor<2, spacedim>> &output) const override;
192 
193  // for documentation, see the Mapping base class
194  virtual void
195  transform(const ArrayView<const Tensor<2, dim>> & input,
196  const MappingType type,
197  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
198  const ArrayView<Tensor<2, spacedim>> &output) const override;
199 
200  // for documentation, see the Mapping base class
201  virtual void
203  const MappingType type,
204  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
205  const ArrayView<Tensor<3, spacedim>> &output) const override;
206 
207  // for documentation, see the Mapping base class
208  virtual void
209  transform(const ArrayView<const Tensor<3, dim>> & input,
210  const MappingType type,
211  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
212  const ArrayView<Tensor<3, spacedim>> &output) const override;
213 
222  unsigned int
223  get_degree() const;
224 
230  get_component_mask() const;
231 
236 
237 private:
243  // documentation can be found in Mapping::requires_update_flags()
244  virtual UpdateFlags
245  requires_update_flags(const UpdateFlags update_flags) const override;
246 
247 public:
259  class InternalData : public Mapping<dim, spacedim>::InternalDataBase
260  {
261  public:
266  const ComponentMask & mask);
267 
272  const double &
273  shape(const unsigned int qpoint, const unsigned int shape_nr) const;
274 
278  double &
279  shape(const unsigned int qpoint, const unsigned int shape_nr);
280 
284  const Tensor<1, dim> &
285  derivative(const unsigned int qpoint, const unsigned int shape_nr) const;
286 
291  derivative(const unsigned int qpoint, const unsigned int shape_nr);
292 
296  const Tensor<2, dim> &
297  second_derivative(const unsigned int qpoint,
298  const unsigned int shape_nr) const;
299 
304  second_derivative(const unsigned int qpoint, const unsigned int shape_nr);
305 
309  const Tensor<3, dim> &
310  third_derivative(const unsigned int qpoint,
311  const unsigned int shape_nr) const;
312 
317  third_derivative(const unsigned int qpoint, const unsigned int shape_nr);
318 
322  const Tensor<4, dim> &
323  fourth_derivative(const unsigned int qpoint,
324  const unsigned int shape_nr) const;
325 
330  fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr);
331 
335  virtual std::size_t
336  memory_consumption() const override;
337 
343  std::vector<double> shape_values;
344 
350  std::vector<Tensor<1, dim>> shape_derivatives;
351 
358  std::vector<Tensor<2, dim>> shape_second_derivatives;
359 
366  std::vector<Tensor<3, dim>> shape_third_derivatives;
367 
374  std::vector<Tensor<4, dim>> shape_fourth_derivatives;
375 
389  std::array<std::vector<Tensor<1, dim>>,
392 
399  unsigned int n_shape_functions;
400 
414 
424  mutable std::vector<DerivativeForm<1, dim, spacedim>> covariant;
425 
433  mutable std::vector<DerivativeForm<1, dim, spacedim>> contravariant;
434 
439  mutable std::vector<double> volume_elements;
440 
444  mutable std::vector<std::vector<Tensor<1, spacedim>>> aux;
445 
449  mutable std::vector<types::global_dof_index> local_dof_indices;
450 
454  mutable std::vector<double> local_dof_values;
455  };
456 
457 private:
458  // documentation can be found in Mapping::get_data()
459  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
460  get_data(const UpdateFlags, const Quadrature<dim> &quadrature) const override;
461 
462  // documentation can be found in Mapping::get_face_data()
463  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
464  get_face_data(const UpdateFlags flags,
465  const Quadrature<dim - 1> &quadrature) const override;
466 
467  // documentation can be found in Mapping::get_subface_data()
468  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
469  get_subface_data(const UpdateFlags flags,
470  const Quadrature<dim - 1> &quadrature) const override;
471 
472  // documentation can be found in Mapping::fill_fe_values()
474  fill_fe_values(
475  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
476  const CellSimilarity::Similarity cell_similarity,
477  const Quadrature<dim> & quadrature,
478  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
480  &output_data) const override;
481 
482  // documentation can be found in Mapping::fill_fe_face_values()
483  virtual void
484  fill_fe_face_values(
485  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
486  const unsigned int face_no,
487  const Quadrature<dim - 1> & quadrature,
488  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
490  &output_data) const override;
491 
492  // documentation can be found in Mapping::fill_fe_subface_values()
493  virtual void
494  fill_fe_subface_values(
495  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
496  const unsigned int face_no,
497  const unsigned int subface_no,
498  const Quadrature<dim - 1> & quadrature,
499  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
501  &output_data) const override;
502 
511  SmartPointer<const VectorType,
514 
518  SmartPointer<const DoFHandlerType,
521 
522 private:
539  do_transform_unit_to_real_cell(const InternalData &mdata) const;
540 
541 
556  Point<dim>
558  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
559  const Point<spacedim> & p,
560  const Point<dim> & initial_p_unit,
561  InternalData & mdata) const;
562 
566  void
568  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
570  InternalData &data) const;
571 
575  virtual void
577  const std::vector<Point<dim>> &unit_points,
579  InternalData &data) const;
580 
581  /*
582  * Which components to use for the mapping.
583  */
584  const ComponentMask fe_mask;
585 
595  std::vector<unsigned int> fe_to_real;
596 
602 
607 
608  void
609  compute_data(const UpdateFlags update_flags,
610  const Quadrature<dim> &q,
611  const unsigned int n_original_q_points,
612  InternalData & data) const;
613 
614  void
615  compute_face_data(const UpdateFlags update_flags,
616  const Quadrature<dim> &q,
617  const unsigned int n_original_q_points,
618  InternalData & data) const;
619 
620 
624  template <int, int, class, class>
625  friend class MappingFEField;
626 };
627 
630 /* -------------- declaration of explicit specializations ------------- */
631 
632 #ifndef DOXYGEN
633 
634 #endif // DOXYGEN
635 
636 DEAL_II_NAMESPACE_CLOSE
637 
638 #endif
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const override
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
SmartPointer< const VectorType, MappingFEField< dim, spacedim, VectorType, DoFHandlerType > > euler_vector
virtual void compute_shapes_virtual(const std::vector< Point< dim >> &unit_points, typename MappingFEField< dim, spacedim, VectorType, DoFHandlerType >::InternalData &data) const
MappingType
Definition: mapping.h:61
std::vector< unsigned int > fe_to_real
Point< spacedim > do_transform_unit_to_real_cell(const InternalData &mdata) const
unsigned int get_degree() const
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
Threads::Mutex fe_values_mutex
std::vector< Tensor< 4, dim > > shape_fourth_derivatives
Abstract base class for mapping classes.
Definition: dof_tools.h:57
virtual void transform(const ArrayView< const Tensor< 1, dim >> &input, const MappingType type, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const ArrayView< Tensor< 1, spacedim >> &output) const override
#define DeclException0(Exception0)
Definition: exceptions.h:473
std::vector< double > volume_elements
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
const Tensor< 2, dim > & second_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
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
Definition: fe.h:38
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const override
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()
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