Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
mapping.cc
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2001 - 2019 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 
17 #include <deal.II/boost_adaptors/bounding_box.h>
18 
19 #include <deal.II/fe/mapping.h>
20 
21 #include <deal.II/grid/tria.h>
22 
23 DEAL_II_NAMESPACE_OPEN
24 
25 
26 template <int dim, int spacedim>
27 std::array<Point<spacedim>, GeometryInfo<dim>::vertices_per_cell>
29  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
30 {
31  std::array<Point<spacedim>, GeometryInfo<dim>::vertices_per_cell> vertices;
32  for (unsigned int i = 0; i < GeometryInfo<dim>::vertices_per_cell; ++i)
33  {
34  vertices[i] = cell->vertex(i);
35  }
36  return vertices;
37 }
38 
39 
40 
41 template <int dim, int spacedim>
45  const bool map_center_of_reference_cell) const
46 {
47  if (map_center_of_reference_cell)
48  {
49  Point<dim> reference_center;
50  for (unsigned int d = 0; d < dim; ++d)
51  reference_center[d] = .5;
52  return transform_unit_to_real_cell(cell, reference_center);
53  }
54  else
55  {
56  const auto vertices = get_vertices(cell);
57  Point<spacedim> center;
58  for (const auto &v : vertices)
59  center += v;
61  }
62 }
63 
64 
65 
66 template <int dim, int spacedim>
69  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const
70 {
71  if (preserves_vertex_locations())
72  return cell->bounding_box();
73  else
74  return BoundingBox<spacedim>(get_vertices(cell));
75 }
76 
77 
78 
79 template <int dim, int spacedim>
80 Point<dim - 1>
83  const unsigned int face_no,
84  const Point<spacedim> & p) const
85 {
86  // The function doesn't make physical sense for dim=1
87  Assert(dim > 1, ExcNotImplemented());
88  // Not implemented for higher dimensions
89  Assert(dim <= 3, ExcNotImplemented());
90 
91  Point<dim> unit_cell_pt = transform_real_to_unit_cell(cell, p);
92 
93  const unsigned int unit_normal_direction =
95 
96  if (dim == 2)
97  {
98  if (unit_normal_direction == 0)
99  return Point<dim - 1>{unit_cell_pt(1)};
100  else if (unit_normal_direction == 1)
101  return Point<dim - 1>{unit_cell_pt(0)};
102  }
103  else if (dim == 3)
104  {
105  if (unit_normal_direction == 0)
106  return Point<dim - 1>{unit_cell_pt(1), unit_cell_pt(2)};
107  else if (unit_normal_direction == 1)
108  return Point<dim - 1>{unit_cell_pt(0), unit_cell_pt(2)};
109  else if (unit_normal_direction == 2)
110  return Point<dim - 1>{unit_cell_pt(0), unit_cell_pt(1)};
111  }
112 
113  // We should never get here
114  Assert(false, ExcInternalError());
115  return {};
116 }
117 
118 /* ---------------------------- InternalDataBase --------------------------- */
119 
120 
121 template <int dim, int spacedim>
123  : update_each(update_default)
124 {}
125 
126 
127 
128 template <int dim, int spacedim>
129 std::size_t
131 {
132  return sizeof(*this);
133 }
134 
135 
136 /*------------------------------ InternalData ------------------------------*/
137 
138 
139 
140 // explicit instantiations
141 #include "mapping.inst"
142 
143 
144 DEAL_II_NAMESPACE_CLOSE
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
Definition: mapping.cc:81
virtual Point< spacedim > get_center(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const bool map_center_of_reference_cell=true) const
Definition: mapping.cc:43
virtual std::size_t memory_consumption() const
Definition: mapping.cc:130
virtual std::array< Point< spacedim >, GeometryInfo< dim >::vertices_per_cell > get_vertices(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const
Definition: mapping.cc:28
No update.
#define Assert(cond, exc)
Definition: exceptions.h:1407
virtual BoundingBox< spacedim > get_bounding_box(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const
Definition: mapping.cc:68
static ::ExceptionBase & ExcNotImplemented()
static ::ExceptionBase & ExcInternalError()