Reference documentation for deal.II version Git ac854ef673 2021-03-06 10:09:45 +0100
\(\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\}}\)
reference_cell.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 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 
20 
24 #include <deal.II/fe/fe_wedge_p.h>
25 #include <deal.II/fe/mapping_fe.h>
26 #include <deal.II/fe/mapping_q1.h>
28 
31 #include <deal.II/grid/tria.h>
32 
33 #include <memory>
34 
36 
37 
38 std::string
40 {
41  if (*this == ReferenceCells::Vertex)
42  return "Vertex";
43  else if (*this == ReferenceCells::Line)
44  return "Line";
45  else if (*this == ReferenceCells::Triangle)
46  return "Tri";
47  else if (*this == ReferenceCells::Quadrilateral)
48  return "Quad";
49  else if (*this == ReferenceCells::Tetrahedron)
50  return "Tet";
51  else if (*this == ReferenceCells::Pyramid)
52  return "Pyramid";
53  else if (*this == ReferenceCells::Wedge)
54  return "Wedge";
55  else if (*this == ReferenceCells::Hexahedron)
56  return "Hex";
57  else if (*this == ReferenceCells::Invalid)
58  return "Invalid";
59 
60  Assert(false, ExcNotImplemented());
61 
62  return "Invalid";
63 }
64 
65 
66 
67 template <int dim, int spacedim>
68 std::unique_ptr<Mapping<dim, spacedim>>
69 ReferenceCell::get_default_mapping(const unsigned int degree) const
70 {
72 
73  if (is_hyper_cube())
74  return std::make_unique<MappingQGeneric<dim, spacedim>>(degree);
75  else if (is_simplex())
76  return std::make_unique<MappingFE<dim, spacedim>>(
78  else if (*this == ReferenceCells::Pyramid)
79  return std::make_unique<MappingFE<dim, spacedim>>(
81  else if (*this == ReferenceCells::Wedge)
82  return std::make_unique<MappingFE<dim, spacedim>>(
83  FE_WedgeP<dim, spacedim>(degree));
84  else
85  {
86  Assert(false, ExcNotImplemented());
87  }
88 
89  return std::make_unique<MappingQGeneric<dim, spacedim>>(degree);
90 }
91 
92 
93 
94 template <int dim, int spacedim>
97 {
99 
100  if (is_hyper_cube())
101  {
103  }
104  else if (is_simplex())
105  {
106  static const MappingFE<dim, spacedim> mapping(
108  return mapping;
109  }
110  else if (*this == ReferenceCells::Pyramid)
111  {
112  static const MappingFE<dim, spacedim> mapping(
114  return mapping;
115  }
116  else if (*this == ReferenceCells::Wedge)
117  {
118  static const MappingFE<dim, spacedim> mapping(
120  return mapping;
121  }
122  else
123  {
124  Assert(false, ExcNotImplemented());
125  }
126 
127  return StaticMappingQ1<dim, spacedim>::mapping; // never reached
128 }
129 
130 
131 
132 template <int dim>
134 ReferenceCell::get_gauss_type_quadrature(const unsigned n_points_1D) const
135 {
137 
138  if (is_hyper_cube())
139  return QGauss<dim>(n_points_1D);
140  else if (is_simplex())
141  return QGaussSimplex<dim>(n_points_1D);
142  else if (*this == ReferenceCells::Pyramid)
143  return QGaussPyramid<dim>(n_points_1D);
144  else if (*this == ReferenceCells::Wedge)
145  return QGaussWedge<dim>(n_points_1D);
146  else
147  Assert(false, ExcNotImplemented());
148 
149  return Quadrature<dim>(); // never reached
150 }
151 
152 
153 
154 template <int dim>
155 const Quadrature<dim> &
157 {
159 
160  // A function that is used to fill a quadrature object of the
161  // desired type the first time we encounter a particular
162  // reference cell
163  const auto create_quadrature = [](const ReferenceCell &reference_cell) {
164  Triangulation<dim> tria;
166 
167  return Quadrature<dim>(tria.get_vertices());
168  };
169 
170  if (is_hyper_cube())
171  {
172  static const Quadrature<dim> quadrature = create_quadrature(*this);
173  return quadrature;
174  }
175  else if (is_simplex())
176  {
177  static const Quadrature<dim> quadrature = create_quadrature(*this);
178  return quadrature;
179  }
180  else if (*this == ReferenceCells::Pyramid)
181  {
182  static const Quadrature<dim> quadrature = create_quadrature(*this);
183  return quadrature;
184  }
185  else if (*this == ReferenceCells::Wedge)
186  {
187  static const Quadrature<dim> quadrature = create_quadrature(*this);
188  return quadrature;
189  }
190  else
191  Assert(false, ExcNotImplemented());
192 
193  static const Quadrature<dim> dummy;
194  return dummy; // never reached
195 }
196 
197 
198 
199 unsigned int
200 ReferenceCell::exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const
201 {
202  AssertIndexRange(vertex_n, n_vertices());
203 
204  if (*this == ReferenceCells::Line)
205  {
206  return vertex_n;
207  }
208  else if (*this == ReferenceCells::Triangle)
209  {
210  return vertex_n;
211  }
212  else if (*this == ReferenceCells::Quadrilateral)
213  {
214  constexpr std::array<unsigned int, 4> exodus_to_deal{{0, 1, 3, 2}};
215  return exodus_to_deal[vertex_n];
216  }
217  else if (*this == ReferenceCells::Tetrahedron)
218  {
219  return vertex_n;
220  }
221  else if (*this == ReferenceCells::Hexahedron)
222  {
223  constexpr std::array<unsigned int, 8> exodus_to_deal{
224  {0, 1, 3, 2, 4, 5, 7, 6}};
225  return exodus_to_deal[vertex_n];
226  }
227  else if (*this == ReferenceCells::Wedge)
228  {
229  constexpr std::array<unsigned int, 6> exodus_to_deal{{2, 1, 0, 5, 4, 3}};
230  return exodus_to_deal[vertex_n];
231  }
232  else if (*this == ReferenceCells::Pyramid)
233  {
234  constexpr std::array<unsigned int, 5> exodus_to_deal{{0, 1, 3, 2, 4}};
235  return exodus_to_deal[vertex_n];
236  }
237 
238  Assert(false, ExcNotImplemented());
239 
241 }
242 
243 
244 
245 unsigned int
246 ReferenceCell::exodusii_face_to_deal_face(const unsigned int face_n) const
247 {
248  AssertIndexRange(face_n, n_faces());
249 
250  if (*this == ReferenceCells::Vertex)
251  {
252  return 0;
253  }
254  if (*this == ReferenceCells::Line)
255  {
256  return face_n;
257  }
258  else if (*this == ReferenceCells::Triangle)
259  {
260  return face_n;
261  }
262  else if (*this == ReferenceCells::Quadrilateral)
263  {
264  constexpr std::array<unsigned int, 4> exodus_to_deal{{2, 1, 3, 0}};
265  return exodus_to_deal[face_n];
266  }
267  else if (*this == ReferenceCells::Tetrahedron)
268  {
269  constexpr std::array<unsigned int, 4> exodus_to_deal{{1, 3, 2, 0}};
270  return exodus_to_deal[face_n];
271  }
272  else if (*this == ReferenceCells::Hexahedron)
273  {
274  constexpr std::array<unsigned int, 6> exodus_to_deal{{2, 1, 3, 0, 4, 5}};
275  return exodus_to_deal[face_n];
276  }
277  else if (*this == ReferenceCells::Wedge)
278  {
279  constexpr std::array<unsigned int, 6> exodus_to_deal{{3, 4, 2, 0, 1}};
280  return exodus_to_deal[face_n];
281  }
282  else if (*this == ReferenceCells::Pyramid)
283  {
284  constexpr std::array<unsigned int, 5> exodus_to_deal{{3, 2, 4, 1, 0}};
285  return exodus_to_deal[face_n];
286  }
287 
288  Assert(false, ExcNotImplemented());
289 
291 }
292 
293 
294 
295 #include "reference_cell.inst"
296 
constexpr const ReferenceCell Pyramid
static const unsigned int invalid_unsigned_int
Definition: types.h:196
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1623
bool is_hyper_cube() const
unsigned int exodusii_vertex_to_deal_vertex(const unsigned int vertex_n) const
unsigned int exodusii_face_to_deal_face(const unsigned int face_n) const
bool is_simplex() const
#define AssertIndexRange(index, range)
Definition: exceptions.h:1691
constexpr const ReferenceCell Triangle
std::string to_string() const
const Quadrature< dim > & get_nodal_type_quadrature() const
constexpr const ReferenceCell Line
constexpr const ReferenceCell Wedge
constexpr const ReferenceCell Tetrahedron
#define Assert(cond, exc)
Definition: exceptions.h:1466
void reference_cell(const ReferenceCell &reference_cell, Triangulation< dim, spacedim > &tria)
Abstract base class for mapping classes.
Definition: mapping.h:303
Quadrature< dim > get_gauss_type_quadrature(const unsigned n_points_1D) const
const std::vector< Point< spacedim > > & get_vertices() const
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:394
constexpr const ReferenceCell Hexahedron
unsigned int n_vertices() const
const Mapping< dim, spacedim > & get_default_linear_mapping() const
unsigned int get_dimension() const
constexpr const ReferenceCell Invalid
std::unique_ptr< Mapping< dim, spacedim > > get_default_mapping(const unsigned int degree) const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:393
constexpr const ReferenceCell Quadrilateral
constexpr const ReferenceCell Vertex
static ::ExceptionBase & ExcNotImplemented()
unsigned int n_faces() const