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
fe_wedge_p.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2020 - 2022 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#include <deal.II/base/config.h>
17
20
21#include <deal.II/fe/fe_dgq.h>
24#include <deal.II/fe/fe_q.h>
27#include <deal.II/fe/fe_tools.h>
29
31
32namespace
33{
38 get_dpo_vector_fe_wedge_p(const unsigned int degree)
39 {
41
42 if (degree == 1)
43 {
44 dpo.dofs_per_object_exclusive = {{1}, {0}, {0, 0, 0, 0, 0}, {0}};
45 dpo.dofs_per_object_inclusive = {{1}, {2}, {3, 3, 4, 4, 4}, {6}};
46 dpo.object_index = {{}, {6}, {6}, {6}};
47 dpo.first_object_index_on_face = {{}, {3, 3, 4, 4, 4}, {3, 3, 4, 4, 4}};
48 }
49 else if (degree == 2)
50 {
51 dpo.dofs_per_object_exclusive = {{1}, {1}, {0, 0, 1, 1, 1}, {0}};
52 dpo.dofs_per_object_inclusive = {{1}, {3}, {6, 6, 9, 9, 9}, {18}};
53 dpo.object_index = {{}, {6}, {15, 15, 15, 16, 17}, {18}};
54 dpo.first_object_index_on_face = {{}, {3, 3, 4, 4, 4}, {6, 6, 8, 8, 8}};
55 }
56 else
57 {
58 Assert(false, ExcNotImplemented());
59 }
60
61 return dpo;
62 }
63
68 get_dpo_vector_fe_wedge_dgp(const unsigned int degree)
69 {
70 unsigned int n_dofs = 0;
71
72 if (degree == 1)
73 n_dofs = 6;
74 else if (degree == 2)
75 n_dofs = 18;
76 else
77 Assert(false, ExcNotImplemented());
78
79 return internal::expand(3, {{0, 0, 0, n_dofs}}, ReferenceCells::Wedge);
80 }
81} // namespace
82
83
84
85template <int dim, int spacedim>
87 const unsigned int degree,
89 const typename FiniteElementData<dim>::Conformity conformity)
90 : ::FE_Poly<dim, spacedim>(
92 FiniteElementData<dim>(dpos,
93 ReferenceCells::Wedge,
94 1,
95 degree,
96 conformity),
97 std::vector<bool>(
98 FiniteElementData<dim>(dpos, ReferenceCells::Wedge, 1, degree)
99 .dofs_per_cell,
100 true),
101 std::vector<ComponentMask>(
102 FiniteElementData<dim>(dpos, ReferenceCells::Wedge, 1, degree)
103 .dofs_per_cell,
104 std::vector<bool>(1, true)))
105{
106 AssertDimension(dim, 3);
107
108 Assert(1 <= degree && degree <= 2, ExcNotImplemented());
109
110 FE_SimplexP<2> fe_triangle(degree);
111 FE_Q<1> fe_line(degree);
112 FE_Q<2> fe_quad(degree);
113
114 for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
115 {
116 const auto pair = this->degree == 1 ? internal::wedge_table_1[i] :
118
119 this->unit_support_points.emplace_back(
120 fe_triangle.get_unit_support_points()[pair[0]][0],
121 fe_triangle.get_unit_support_points()[pair[0]][1],
122 fe_line.get_unit_support_points()[pair[1]][0]);
123 }
124
125 this->unit_face_support_points.resize(this->reference_cell().n_faces());
126
127 for (const auto f : this->reference_cell().face_indices())
128 if (this->reference_cell().face_reference_cell(f) ==
130 for (const auto &p : fe_triangle.get_unit_support_points())
131 this->unit_face_support_points[f].emplace_back(p[0], p[1]);
132 else if (this->reference_cell().face_reference_cell(f) ==
134 for (const auto &p : fe_quad.get_unit_support_points())
135 this->unit_face_support_points[f].emplace_back(p[0], p[1]);
136 else
137 Assert(false, ExcInternalError());
138}
139
140
141
142template <int dim, int spacedim>
143void
146 const std::vector<Vector<double>> &support_point_values,
147 std::vector<double> & nodal_values) const
148{
149 AssertDimension(support_point_values.size(),
150 this->get_unit_support_points().size());
151 AssertDimension(support_point_values.size(), nodal_values.size());
152 AssertDimension(this->dofs_per_cell, nodal_values.size());
153
154 for (unsigned int i = 0; i < this->dofs_per_cell; ++i)
155 {
156 AssertDimension(support_point_values[i].size(), 1);
157
158 nodal_values[i] = support_point_values[i](0);
159 }
160}
161
162
163
164template <int dim, int spacedim>
165FE_WedgeP<dim, spacedim>::FE_WedgeP(const unsigned int degree)
166 : FE_WedgePoly<dim, spacedim>(degree,
167 get_dpo_vector_fe_wedge_p(degree),
168 FiniteElementData<dim>::H1)
169{}
170
171
172
173template <int dim, int spacedim>
174std::unique_ptr<FiniteElement<dim, spacedim>>
176{
177 return std::make_unique<FE_WedgeP<dim, spacedim>>(*this);
178}
179
180
181
182template <int dim, int spacedim>
183std::string
185{
186 std::ostringstream namebuf;
187 namebuf << "FE_WedgeP<" << dim << ">(" << this->degree << ")";
188
189 return namebuf.str();
190}
191
192
193
194template <int dim, int spacedim>
197 const FiniteElement<dim, spacedim> &fe_other,
198 const unsigned int codim) const
199{
200 Assert(codim <= dim, ExcImpossibleInDim(dim));
201
202 // vertex/line/face domination
203 // (if fe_other is derived from FE_SimplexDGP)
204 // ------------------------------------
205 if (codim > 0)
206 if (dynamic_cast<const FE_SimplexDGP<dim, spacedim> *>(&fe_other) !=
207 nullptr)
208 // there are no requirements between continuous and discontinuous
209 // elements
211
212
213 // vertex/line/face domination
214 // (if fe_other is not derived from FE_SimplexDGP)
215 // & cell domination
216 // ----------------------------------------
217 if (const FE_WedgeP<dim, spacedim> *fe_wp_other =
218 dynamic_cast<const FE_WedgeP<dim, spacedim> *>(&fe_other))
219 {
220 if (this->degree < fe_wp_other->degree)
222 else if (this->degree == fe_wp_other->degree)
224 else
226 }
227 else if (const FE_SimplexP<dim, spacedim> *fe_p_other =
228 dynamic_cast<const FE_SimplexP<dim, spacedim> *>(&fe_other))
229 {
230 if (this->degree < fe_p_other->degree)
232 else if (this->degree == fe_p_other->degree)
234 else
236 }
237 else if (const FE_Q<dim, spacedim> *fe_q_other =
238 dynamic_cast<const FE_Q<dim, spacedim> *>(&fe_other))
239 {
240 if (this->degree < fe_q_other->degree)
242 else if (this->degree == fe_q_other->degree)
244 else
246 }
247 else if (const FE_Nothing<dim> *fe_nothing =
248 dynamic_cast<const FE_Nothing<dim> *>(&fe_other))
249 {
250 if (fe_nothing->is_dominating())
252 else
253 // the FE_Nothing has no degrees of freedom and it is typically used
254 // in a context where we don't require any continuity along the
255 // interface
257 }
258
259 Assert(false, ExcNotImplemented());
261}
262
263
264
265template <int dim, int spacedim>
266std::vector<std::pair<unsigned int, unsigned int>>
268 const FiniteElement<dim, spacedim> &fe_other) const
269{
270 (void)fe_other;
271
272 Assert((dynamic_cast<const FE_SimplexP<dim, spacedim> *>(&fe_other)) ||
273 (dynamic_cast<const FE_Q<dim, spacedim> *>(&fe_other)),
275
276 return {{0, 0}};
277}
278
279
280
281template <int dim, int spacedim>
282std::vector<std::pair<unsigned int, unsigned int>>
284 const FiniteElement<dim, spacedim> &fe_other) const
285{
286 (void)fe_other;
287
288 Assert((dynamic_cast<const FE_SimplexP<dim, spacedim> *>(&fe_other)) ||
289 (dynamic_cast<const FE_Q<dim, spacedim> *>(&fe_other)),
291
292 std::vector<std::pair<unsigned int, unsigned int>> result;
293
294 for (unsigned int i = 0; i < this->degree - 1; ++i)
295 result.emplace_back(i, i);
296
297 return result;
298}
299
300
301
302template <int dim, int spacedim>
303std::vector<std::pair<unsigned int, unsigned int>>
305 const FiniteElement<dim, spacedim> &fe_other,
306 const unsigned int face_no) const
307{
308 (void)fe_other;
309
310 AssertIndexRange(face_no, 5);
311
312 if (face_no < 2)
313 {
314 Assert((dynamic_cast<const FE_SimplexP<dim, spacedim> *>(&fe_other)),
316 }
317 else
318 {
319 Assert((dynamic_cast<const FE_Q<dim, spacedim> *>(&fe_other)),
321 }
322
323 std::vector<std::pair<unsigned int, unsigned int>> result;
324
325 for (unsigned int i = 0; i < this->n_dofs_per_quad(face_no); ++i)
326 result.emplace_back(i, i);
327
328 return result;
329}
330
331
332
333template <int dim, int spacedim>
335 : FE_WedgePoly<dim, spacedim>(degree,
336 get_dpo_vector_fe_wedge_dgp(degree),
337 FiniteElementData<dim>::L2)
338{}
339
340
341
342template <int dim, int spacedim>
343std::unique_ptr<FiniteElement<dim, spacedim>>
345{
346 return std::make_unique<FE_WedgeDGP<dim, spacedim>>(*this);
347}
348
349
350
351template <int dim, int spacedim>
352std::string
354{
355 std::ostringstream namebuf;
356 namebuf << "FE_WedgeDGP<" << dim << ">(" << this->degree << ")";
357
358 return namebuf.str();
359}
360
361// explicit instantiations
362#include "fe_wedge_p.inst"
363
Definition fe_q.h:551
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:86
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:453
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:2447
const std::vector< Point< dim > > & get_unit_support_points() const
std::vector< Point< dim > > unit_support_points
Definition fe.h:2440
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:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
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:25
STL namespace.
std::vector< std::vector< unsigned int > > object_index
Definition fe_data.h:198
std::vector< std::vector< unsigned int > > first_object_index_on_face
Definition fe_data.h:203
std::vector< std::vector< unsigned int > > dofs_per_object_inclusive
Definition fe_data.h:193
std::vector< std::vector< unsigned int > > dofs_per_object_exclusive
Definition fe_data.h:188