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.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2001 - 2023 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
18
19#include <deal.II/fe/mapping.h>
20
21#include <deal.II/grid/tria.h>
22
23#ifdef DEAL_II_BOOST_HAS_BROKEN_HEADER_DEPRECATIONS
24# define BOOST_ALLOW_DEPRECATED_HEADERS
25#endif
26#include <boost/geometry.hpp>
27#ifdef DEAL_II_BOOST_HAS_BROKEN_HEADER_DEPRECATIONS
28# undef BOOST_ALLOW_DEPRECATED_HEADERS
29#endif
30
31#include <limits>
32
34
35
36template <int dim, int spacedim>
37boost::container::small_vector<Point<spacedim>,
40 const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
41{
42 boost::container::small_vector<Point<spacedim>,
45 for (const unsigned int i : cell->vertex_indices())
46 vertices.push_back(cell->vertex(i));
47
48 return vertices;
49}
50
51
52
53template <int dim, int spacedim>
54boost::container::small_vector<Point<spacedim>,
58 const unsigned int face_no) const
59{
60 boost::container::small_vector<Point<spacedim>,
62 face_vertices;
63
64 const auto &cell_vertices = get_vertices(cell);
65 const auto &reference_cell = cell->reference_cell();
66 const auto face_orientation = cell->combined_face_orientation(face_no);
67
68 for (const unsigned int v :
69 reference_cell.face_reference_cell(face_no).vertex_indices())
70 {
71 face_vertices.push_back(
72 cell_vertices[reference_cell.face_to_cell_vertices(
73 face_no, v, face_orientation)]);
74 }
75
76 return face_vertices;
77}
78
79
80
81template <int dim, int spacedim>
85 const bool map_barycenter_of_reference_cell) const
86{
87 if (map_barycenter_of_reference_cell)
88 {
89 return transform_unit_to_real_cell(
90 cell, cell->reference_cell().template barycenter<dim>());
91 }
92 else
93 {
94 const auto vertices = get_vertices(cell);
96 for (const auto &v : vertices)
97 center += v;
98 return center / cell->n_vertices();
99 }
100}
101
102
103
104template <int dim, int spacedim>
107 const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
108{
109 if (preserves_vertex_locations())
110 return cell->bounding_box();
111 else
112 return BoundingBox<spacedim>(get_vertices(cell));
113}
114
115
116
117template <int dim, int spacedim>
118void
124{
126}
127
128
129
130template <int dim, int spacedim>
131void
134 const ArrayView<const Point<spacedim>> & real_points,
135 const ArrayView<Point<dim>> & unit_points) const
136{
137 AssertDimension(real_points.size(), unit_points.size());
138 for (unsigned int i = 0; i < real_points.size(); ++i)
139 {
140 try
141 {
142 unit_points[i] = transform_real_to_unit_cell(cell, real_points[i]);
143 }
145 {
146 unit_points[i] = Point<dim>();
147 unit_points[i][0] = std::numeric_limits<double>::infinity();
148 }
149 }
150}
151
152
153
154template <int dim, int spacedim>
155Point<dim - 1>
158 const unsigned int face_no,
159 const Point<spacedim> & p) const
160{
161 // The function doesn't make physical sense for dim=1
162 Assert(dim > 1, ExcNotImplemented());
163 // Not implemented for higher dimensions
164 Assert(dim <= 3, ExcNotImplemented());
165
166 Point<dim> unit_cell_pt = transform_real_to_unit_cell(cell, p);
167
168 const unsigned int unit_normal_direction =
170
171 if (dim == 2)
172 {
173 if (unit_normal_direction == 0)
174 return Point<dim - 1>{unit_cell_pt(1)};
175 else if (unit_normal_direction == 1)
176 return Point<dim - 1>{unit_cell_pt(0)};
177 }
178 else if (dim == 3)
179 {
180 if (unit_normal_direction == 0)
181 return Point<dim - 1>{unit_cell_pt(1), unit_cell_pt(2)};
182 else if (unit_normal_direction == 1)
183 return Point<dim - 1>{unit_cell_pt(0), unit_cell_pt(2)};
184 else if (unit_normal_direction == 2)
185 return Point<dim - 1>{unit_cell_pt(0), unit_cell_pt(1)};
186 }
187
188 // We should never get here
189 Assert(false, ExcInternalError());
190 return {};
191}
192
193
194
195#ifndef DOXYGEN
196template <int dim, int spacedim>
197void
200 const unsigned int face_no,
201 const hp::QCollection<dim - 1> & quadrature,
202 const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
204 &output_data) const
205{
206 // base class version, implement overridden function in derived classes
207 AssertDimension(quadrature.size(), 1);
208 fill_fe_face_values(cell, face_no, quadrature[0], internal_data, output_data);
209}
210
211
212
213template <int dim, int spacedim>
214void
217 const unsigned int face_no,
218 const Quadrature<dim - 1> & quadrature,
219 const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
221 &output_data) const
222{
223 Assert(false,
224 ExcMessage("Use of a deprecated interface, please implement "
225 "fill_fe_face_values taking a hp::QCollection argument"));
226 (void)cell;
227 (void)face_no;
228 (void)quadrature;
229 (void)internal_data;
230 (void)output_data;
231}
232
233
234
235template <int dim, int spacedim>
236std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
238 const UpdateFlags update_flags,
239 const hp::QCollection<dim - 1> &quadrature) const
240{
241 // base class version, implement overridden function in derived classes
242 return get_face_data(update_flags, quadrature[0]);
243}
244
245
246
247template <int dim, int spacedim>
248std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
250 const UpdateFlags update_flags,
251 const Quadrature<dim - 1> &quadrature) const
252{
253 Assert(false,
254 ExcMessage("Use of a deprecated interface, please implement "
255 "fill_fe_face_values taking a hp::QCollection argument"));
256 (void)update_flags;
257 (void)quadrature;
258
259 return std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>();
260}
261#endif
262
263/* ---------------------------- InternalDataBase --------------------------- */
264
265
266template <int dim, int spacedim>
268 : update_each(update_default)
269{}
270
271
272
273template <int dim, int spacedim>
274std::size_t
276{
277 return sizeof(*this);
278}
279
280
281/* ------------------------------ Global functions ------------------------- */
282
283template <int dim, int spacedim>
286{
287 const auto &reference_cells = triangulation.get_reference_cells();
288 Assert(reference_cells.size() == 1,
290 "This function can only work for triangulations that "
291 "use only a single cell type -- for example, only triangles "
292 "or only quadrilaterals. For mixed meshes, there is no "
293 "single linear mapping object that can be used for all "
294 "cells of the triangulation. The triangulation you are "
295 "passing to this function uses multiple cell types."));
296
297 return reference_cells.front()
298 .template get_default_linear_mapping<dim, spacedim>();
299}
300
301
302
303// explicit instantiations
304#include "mapping.inst"
305
306
virtual std::size_t memory_consumption() const
Abstract base class for mapping classes.
Definition mapping.h:317
Point< dim - 1 > project_real_point_to_unit_point_on_face(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const Point< spacedim > &p) const
virtual BoundingBox< spacedim > get_bounding_box(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const
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
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
virtual void transform_points_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const ArrayView< const Point< spacedim > > &real_points, const ArrayView< Point< dim > > &unit_points) const
virtual boost::container::small_vector< Point< spacedim >, GeometryInfo< dim >::vertices_per_cell > get_vertices(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const
virtual Point< spacedim > get_center(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const bool map_barycenter_of_reference_cell=true) const
virtual std::unique_ptr< InternalDataBase > get_face_data(const UpdateFlags update_flags, const hp::QCollection< dim - 1 > &quadrature) const
Definition point.h:112
unsigned int size() const
Definition collection.h:265
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
Point< 3 > center
Point< 3 > vertices[4]
unsigned int vertex_indices[2]
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
UpdateFlags
@ update_default
No update.
const Mapping< dim, spacedim > & get_default_linear_mapping(const Triangulation< dim, spacedim > &triangulation)
Definition mapping.cc:285
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
Point< 1 > transform_real_to_unit_cell(const std::array< Point< spacedim >, GeometryInfo< 1 >::vertices_per_cell > &vertices, const Point< spacedim > &p)
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation