Reference documentation for deal.II version Git 24041777c5 2021-12-02 20:43:25 -0700
\(\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_cartesian.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2001 - 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_cartesian_h
17 #define dealii_mapping_cartesian_h
18 
19 
20 #include <deal.II/base/config.h>
21 
23 
24 #include <deal.II/fe/mapping.h>
25 
26 #include <cmath>
27 
28 
30 
33 
76 template <int dim, int spacedim = dim>
77 class MappingCartesian : public Mapping<dim, spacedim>
78 {
79 public:
80  // for documentation, see the Mapping base class
81  virtual std::unique_ptr<Mapping<dim, spacedim>>
82  clone() const override;
83 
88  virtual bool
89  preserves_vertex_locations() const override;
90 
91  virtual bool
92  is_compatible_with(const ReferenceCell &reference_cell) const override;
93 
99  // for documentation, see the Mapping base class
100  virtual Point<spacedim>
102  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
103  const Point<dim> &p) const override;
104 
105  // for documentation, see the Mapping base class
106  virtual Point<dim>
108  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
109  const Point<spacedim> &p) const override;
110 
120  // for documentation, see the Mapping base class
121  virtual void
122  transform(const ArrayView<const Tensor<1, dim>> & input,
123  const MappingKind kind,
124  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
125  const ArrayView<Tensor<1, spacedim>> &output) const override;
126 
127  // for documentation, see the Mapping base class
128  virtual void
130  const MappingKind kind,
131  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
132  const ArrayView<Tensor<2, spacedim>> &output) const override;
133 
134  // for documentation, see the Mapping base class
135  virtual void
136  transform(const ArrayView<const Tensor<2, dim>> & input,
137  const MappingKind kind,
138  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
139  const ArrayView<Tensor<2, spacedim>> &output) const override;
140 
141  // for documentation, see the Mapping base class
142  virtual void
144  const MappingKind kind,
145  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
146  const ArrayView<Tensor<3, spacedim>> &output) const override;
147 
148  // for documentation, see the Mapping base class
149  virtual void
150  transform(const ArrayView<const Tensor<3, dim>> & input,
151  const MappingKind kind,
152  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
153  const ArrayView<Tensor<3, spacedim>> &output) const override;
154 
176  void
178  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
179  const ArrayView<const Point<dim>> & unit_points,
180  const UpdateFlags update_flags,
182  &output_data) const;
183 
200  class InternalData : public Mapping<dim, spacedim>::InternalDataBase
201  {
202  public:
206  InternalData() = default;
207 
211  InternalData(const Quadrature<dim> &quadrature);
212 
216  virtual std::size_t
217  memory_consumption() const override;
218 
224 
228  mutable double volume_element;
229 
233  std::vector<Point<dim>> quadrature_points;
234  };
235 
236 private:
237  // documentation can be found in Mapping::requires_update_flags()
238  virtual UpdateFlags
239  requires_update_flags(const UpdateFlags update_flags) const override;
240 
241  // documentation can be found in Mapping::get_data()
242  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
243  get_data(const UpdateFlags, const Quadrature<dim> &quadrature) const override;
244 
246 
247  // documentation can be found in Mapping::get_face_data()
248  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
249  get_face_data(const UpdateFlags flags,
250  const hp::QCollection<dim - 1> &quadrature) const override;
251 
252  // documentation can be found in Mapping::get_subface_data()
253  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
254  get_subface_data(const UpdateFlags flags,
255  const Quadrature<dim - 1> &quadrature) const override;
256 
257  // documentation can be found in Mapping::fill_fe_values()
260  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
261  const CellSimilarity::Similarity cell_similarity,
262  const Quadrature<dim> & quadrature,
263  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
265  &output_data) const override;
266 
268 
269  // documentation can be found in Mapping::fill_fe_face_values()
270  virtual void
272  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
273  const unsigned int face_no,
274  const hp::QCollection<dim - 1> & quadrature,
275  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
277  &output_data) const override;
278 
279  // documentation can be found in Mapping::fill_fe_subface_values()
280  virtual void
282  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
283  const unsigned int face_no,
284  const unsigned int subface_no,
285  const Quadrature<dim - 1> & quadrature,
286  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
288  &output_data) const override;
289 
290  // documentation can be found in Mapping::fill_fe_immersed_surface_values()
291  virtual void
293  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
295  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
297  &output_data) const override;
298 
307  void
309  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
310  const CellSimilarity::Similarity cell_similarity,
311  const InternalData & data) const;
312 
319  void
321  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
322  const InternalData & data,
323  std::vector<Point<dim>> &quadrature_points) const;
324 
331  void
333  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
334  const unsigned int face_no,
335  const InternalData & data,
336  std::vector<Point<dim>> &quadrature_points) const;
337 
344  void
346  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
347  const unsigned int face_no,
348  const unsigned int sub_no,
349  const InternalData & data,
350  std::vector<Point<dim>> &quadrature_points) const;
351 
358  void
360  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
361  const InternalData & data,
362  const typename QProjector<dim>::DataSetDescriptor & offset,
363  std::vector<Point<dim>> &quadrature_points) const;
364 
369  void
371  const unsigned int face_no,
372  const InternalData & data,
373  std::vector<Tensor<1, dim>> &normal_vectors) const;
374 
380  void
382  const InternalData & data,
383  const CellSimilarity::Similarity cell_similarity,
385  &output_data) const;
386 
387 
392  void
393  maybe_update_volume_elements(const InternalData &data) const;
394 
399  void
401  const InternalData & data,
402  const CellSimilarity::Similarity cell_similarity,
404  &output_data) const;
405 
410  void
412  const InternalData & data,
413  const CellSimilarity::Similarity cell_similarity,
415  &output_data) const;
416 };
417 
421 
422 #endif
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const override
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 maybe_update_normal_vectors(const unsigned int face_no, const InternalData &data, std::vector< Tensor< 1, dim >> &normal_vectors) const
void maybe_update_volume_elements(const InternalData &data) const
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 update_cell_extents(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const InternalData &data) const
void transform_quadrature_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const InternalData &data, const typename QProjector< dim >::DataSetDescriptor &offset, std::vector< Point< dim >> &quadrature_points) 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
std::vector< Point< dim > > quadrature_points
void maybe_update_face_quadrature_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const InternalData &data, std::vector< Point< dim >> &quadrature_points) const
MappingKind
Definition: mapping.h:71
void fill_mapping_data_for_generic_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const ArrayView< const Point< dim >> &unit_points, const UpdateFlags update_flags, ::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const
virtual std::size_t memory_consumption() const override
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
UpdateFlags
Abstract base class for mapping classes.
Definition: mapping.h:310
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags flags, const Quadrature< dim - 1 > &quadrature) const override
void maybe_update_jacobian_derivatives(const InternalData &data, const CellSimilarity::Similarity cell_similarity, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const
void maybe_update_inverse_jacobians(const InternalData &data, const CellSimilarity::Similarity cell_similarity, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:404
void maybe_update_jacobians(const InternalData &data, const CellSimilarity::Similarity cell_similarity, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const
void maybe_update_subface_quadrature_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const InternalData &data, std::vector< Point< dim >> &quadrature_points) const
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
void maybe_update_cell_quadrature_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const InternalData &data, std::vector< Point< dim >> &quadrature_points) const
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const override
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:403
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
virtual bool preserves_vertex_locations() const override
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const override
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
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
virtual bool is_compatible_with(const ReferenceCell &reference_cell) 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