Reference documentation for deal.II version GIT relicensing-216-gcda843ca70 2024-03-28 12:20:02+00:00
\(\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
fe_wedge_p.cc
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2021 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#include <deal.II/base/config.h>
16
19
20#include <deal.II/fe/fe_dgq.h>
23#include <deal.II/fe/fe_q.h>
26#include <deal.II/fe/fe_tools.h>
28
30
31namespace
32{
37 get_dpo_vector_fe_wedge_p(const unsigned int degree)
38 {
40
41 if (degree == 1)
42 {
43 dpo.dofs_per_object_exclusive = {{1}, {0}, {0, 0, 0, 0, 0}, {0}};
44 dpo.dofs_per_object_inclusive = {{1}, {2}, {3, 3, 4, 4, 4}, {6}};
45 dpo.object_index = {{}, {6}, {6}, {6}};
46 dpo.first_object_index_on_face = {{}, {3, 3, 4, 4, 4}, {3, 3, 4, 4, 4}};
47 }
48 else if (degree == 2)
49 {
50 dpo.dofs_per_object_exclusive = {{1}, {1}, {0, 0, 1, 1, 1}, {0}};
51 dpo.dofs_per_object_inclusive = {{1}, {3}, {6, 6, 9, 9, 9}, {18}};
52 dpo.object_index = {{}, {6}, {15, 15, 15, 16, 17}, {18}};
53 dpo.first_object_index_on_face = {{}, {3, 3, 4, 4, 4}, {6, 6, 8, 8, 8}};
54 }
55 else
56 {
58 }
59
60 return dpo;
61 }
62
67 get_dpo_vector_fe_wedge_dgp(const unsigned int degree)
68 {
69 unsigned int n_dofs = 0;
70
71 if (degree == 1)
72 n_dofs = 6;
73 else if (degree == 2)
74 n_dofs = 18;
75 else
77
78 return internal::expand(3, {{0, 0, 0, n_dofs}}, ReferenceCells::Wedge);
79 }
80} // namespace
81
82
83
84template <int dim, int spacedim>
86 const unsigned int degree,
88 const typename FiniteElementData<dim>::Conformity conformity)
89 : ::FE_Poly<dim, spacedim>(
91 FiniteElementData<dim>(dpos,
92 ReferenceCells::Wedge,
93 1,
94 degree,
95 conformity),
96 std::vector<bool>(
97 FiniteElementData<dim>(dpos, ReferenceCells::Wedge, 1, degree)
98 .dofs_per_cell,
99 true),
100 std::vector<ComponentMask>(
101 FiniteElementData<dim>(dpos, ReferenceCells::Wedge, 1, degree)
102 .dofs_per_cell,
103 ComponentMask(std::vector<bool>(1, true))))
104{
105 AssertDimension(dim, 3);
106
107 Assert(1 <= degree && degree <= 2, ExcNotImplemented());
108
109 FE_SimplexP<2> fe_triangle(degree);
110 FE_Q<1> fe_line(degree);
111 FE_Q<2> fe_quad(degree);
112
113 for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
114 {
115 const auto pair = this->degree == 1 ? internal::wedge_table_1[i] :
117
118 this->unit_support_points.emplace_back(
119 fe_triangle.get_unit_support_points()[pair[0]][0],
120 fe_triangle.get_unit_support_points()[pair[0]][1],
121 fe_line.get_unit_support_points()[pair[1]][0]);
122 }
123
124 this->unit_face_support_points.resize(this->reference_cell().n_faces());
125
126 for (const auto f : this->reference_cell().face_indices())
127 if (this->reference_cell().face_reference_cell(f) ==
129 for (const auto &p : fe_triangle.get_unit_support_points())
130 this->unit_face_support_points[f].emplace_back(p[0], p[1]);
131 else if (this->reference_cell().face_reference_cell(f) ==
133 for (const auto &p : fe_quad.get_unit_support_points())
134 this->unit_face_support_points[f].emplace_back(p[0], p[1]);
135 else
137}
138
139
140
141template <int dim, int spacedim>
142void
145 const std::vector<Vector<double>> &support_point_values,
146 std::vector<double> &nodal_values) const
147{
148 AssertDimension(support_point_values.size(),
149 this->get_unit_support_points().size());
150 AssertDimension(support_point_values.size(), nodal_values.size());
151 AssertDimension(this->dofs_per_cell, nodal_values.size());
152
153 for (unsigned int i = 0; i < this->dofs_per_cell; ++i)
154 {
155 AssertDimension(support_point_values[i].size(), 1);
156
157 nodal_values[i] = support_point_values[i](0);
158 }
159}
160
161
162
163template <int dim, int spacedim>
164FE_WedgeP<dim, spacedim>::FE_WedgeP(const unsigned int degree)
165 : FE_WedgePoly<dim, spacedim>(degree,
166 get_dpo_vector_fe_wedge_p(degree),
167 FiniteElementData<dim>::H1)
168{}
169
170
171
172template <int dim, int spacedim>
173std::unique_ptr<FiniteElement<dim, spacedim>>
175{
176 return std::make_unique<FE_WedgeP<dim, spacedim>>(*this);
177}
178
179
180
181template <int dim, int spacedim>
182std::string
184{
185 std::ostringstream namebuf;
186 namebuf << "FE_WedgeP<" << dim << ">(" << this->degree << ")";
187
188 return namebuf.str();
189}
190
191
192
193template <int dim, int spacedim>
196 const FiniteElement<dim, spacedim> &fe_other,
197 const unsigned int codim) const
198{
199 Assert(codim <= dim, ExcImpossibleInDim(dim));
200
201 // vertex/line/face domination
202 // (if fe_other is derived from FE_SimplexDGP)
203 // ------------------------------------
204 if (codim > 0)
205 if (dynamic_cast<const FE_SimplexDGP<dim, spacedim> *>(&fe_other) !=
206 nullptr)
207 // there are no requirements between continuous and discontinuous
208 // elements
210
211
212 // vertex/line/face domination
213 // (if fe_other is not derived from FE_SimplexDGP)
214 // & cell domination
215 // ----------------------------------------
216 if (const FE_WedgeP<dim, spacedim> *fe_wp_other =
217 dynamic_cast<const FE_WedgeP<dim, spacedim> *>(&fe_other))
218 {
219 if (this->degree < fe_wp_other->degree)
221 else if (this->degree == fe_wp_other->degree)
223 else
225 }
226 else if (const FE_SimplexP<dim, spacedim> *fe_p_other =
227 dynamic_cast<const FE_SimplexP<dim, spacedim> *>(&fe_other))
228 {
229 if (this->degree < fe_p_other->degree)
231 else if (this->degree == fe_p_other->degree)
233 else
235 }
236 else if (const FE_Q<dim, spacedim> *fe_q_other =
237 dynamic_cast<const FE_Q<dim, spacedim> *>(&fe_other))
238 {
239 if (this->degree < fe_q_other->degree)
241 else if (this->degree == fe_q_other->degree)
243 else
245 }
246 else if (const FE_Nothing<dim> *fe_nothing =
247 dynamic_cast<const FE_Nothing<dim> *>(&fe_other))
248 {
249 if (fe_nothing->is_dominating())
251 else
252 // the FE_Nothing has no degrees of freedom and it is typically used
253 // in a context where we don't require any continuity along the
254 // interface
256 }
257
260}
261
262
263
264template <int dim, int spacedim>
265std::vector<std::pair<unsigned int, unsigned int>>
267 const FiniteElement<dim, spacedim> &fe_other) const
268{
269 (void)fe_other;
270
271 Assert((dynamic_cast<const FE_SimplexP<dim, spacedim> *>(&fe_other)) ||
272 (dynamic_cast<const FE_Q<dim, spacedim> *>(&fe_other)),
274
275 return {{0, 0}};
276}
277
278
279
280template <int dim, int spacedim>
281std::vector<std::pair<unsigned int, unsigned int>>
283 const FiniteElement<dim, spacedim> &fe_other) const
284{
285 (void)fe_other;
286
287 Assert((dynamic_cast<const FE_SimplexP<dim, spacedim> *>(&fe_other)) ||
288 (dynamic_cast<const FE_Q<dim, spacedim> *>(&fe_other)),
290
291 std::vector<std::pair<unsigned int, unsigned int>> result;
292
293 for (unsigned int i = 0; i < this->degree - 1; ++i)
294 result.emplace_back(i, i);
295
296 return result;
297}
298
299
300
301template <int dim, int spacedim>
302std::vector<std::pair<unsigned int, unsigned int>>
304 const FiniteElement<dim, spacedim> &fe_other,
305 const unsigned int face_no) const
306{
307 (void)fe_other;
308
309 AssertIndexRange(face_no, 5);
310
311 if (face_no < 2)
312 {
313 Assert((dynamic_cast<const FE_SimplexP<dim, spacedim> *>(&fe_other)),
315 }
316 else
317 {
318 Assert((dynamic_cast<const FE_Q<dim, spacedim> *>(&fe_other)),
320 }
321
322 std::vector<std::pair<unsigned int, unsigned int>> result;
323
324 for (unsigned int i = 0; i < this->n_dofs_per_quad(face_no); ++i)
325 result.emplace_back(i, i);
326
327 return result;
328}
329
330
331
332template <int dim, int spacedim>
334 : FE_WedgePoly<dim, spacedim>(degree,
335 get_dpo_vector_fe_wedge_dgp(degree),
336 FiniteElementData<dim>::L2)
337{}
338
339
340
341template <int dim, int spacedim>
342std::unique_ptr<FiniteElement<dim, spacedim>>
344{
345 return std::make_unique<FE_WedgeDGP<dim, spacedim>>(*this);
346}
347
348
349
350template <int dim, int spacedim>
351std::string
353{
354 std::ostringstream namebuf;
355 namebuf << "FE_WedgeDGP<" << dim << ">(" << this->degree << ")";
356
357 return namebuf.str();
358}
359
360// explicit instantiations
361#include "fe_wedge_p.inst"
362
Definition fe_q.h:550
std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
std::string get_name() const override
FE_WedgeDGP(const unsigned int degree)
std::vector< std::pair< unsigned int, unsigned int > > hp_quad_dof_identities(const FiniteElement< dim, spacedim > &fe_other, const unsigned int face_no=0) const override
FE_WedgeP(const unsigned int degree)
std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim) const override
std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
std::string get_name() const override
FE_WedgePoly(const unsigned int degree, const internal::GenericDoFsPerObject &dpos, const typename FiniteElementData< dim >::Conformity conformity)
Definition fe_wedge_p.cc:85
virtual void convert_generalized_support_point_values_to_dof_values(const std::vector< Vector< double > > &support_point_values, std::vector< double > &nodal_values) const override
const unsigned int degree
Definition fe_data.h:452
unsigned int n_dofs_per_cell() const
ReferenceCell reference_cell() const
std::vector< std::vector< Point< dim - 1 > > > unit_face_support_points
Definition fe.h:2442
const std::vector< Point< dim > > & get_unit_support_points() const
std::vector< Point< dim > > unit_support_points
Definition fe.h:2435
unsigned int n_faces() const
ReferenceCell face_reference_cell(const unsigned int face_no) const
std_cxx20::ranges::iota_view< unsigned int, unsigned int > face_indices() const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
#define DEAL_II_ASSERT_UNREACHABLE()
#define DEAL_II_NOT_IMPLEMENTED()
constexpr const ReferenceCell Quadrilateral
constexpr const ReferenceCell Wedge
constexpr const ReferenceCell Triangle
static const constexpr ::ndarray< unsigned int, 18, 2 > wedge_table_2
static const constexpr ::ndarray< unsigned int, 6, 2 > wedge_table_1
internal::GenericDoFsPerObject expand(const unsigned int dim, const std::vector< unsigned int > &dofs_per_object, const ReferenceCell reference_cell)
Definition fe_data.cc:24
STL namespace.
std::vector< std::vector< unsigned int > > object_index
Definition fe_data.h:197
std::vector< std::vector< unsigned int > > first_object_index_on_face
Definition fe_data.h:202
std::vector< std::vector< unsigned int > > dofs_per_object_inclusive
Definition fe_data.h:192
std::vector< std::vector< unsigned int > > dofs_per_object_exclusive
Definition fe_data.h:187