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_trace.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2000 - 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#include <deal.II/base/config.h>
17
20
23#include <deal.II/fe/fe_q.h>
24#include <deal.II/fe/fe_tools.h>
25#include <deal.II/fe/fe_trace.h>
26
27#include <fstream>
28#include <iostream>
29#include <memory>
30#include <sstream>
31
33
34
35
36template <int dim, int spacedim>
37FE_TraceQ<dim, spacedim>::FE_TraceQ(const unsigned int degree)
38 : FE_PolyFace<TensorProductPolynomials<dim - 1>, dim, spacedim>(
40 Polynomials::generate_complete_Lagrange_basis(
41 QGaussLobatto<1>(degree + 1).get_points())),
42 FiniteElementData<dim>(get_dpo_vector(degree),
43 1,
44 degree,
45 FiniteElementData<dim>::L2),
46 std::vector<bool>(1, true))
47 , fe_q(degree)
48{
49 Assert(degree > 0,
50 ExcMessage("FE_Trace can only be used for polynomial degrees "
51 "greater than zero"));
52 this->poly_space.set_numbering(
53 FETools::hierarchic_to_lexicographic_numbering<dim - 1>(degree));
54
55 // Initialize face support points
56 AssertDimension(this->n_unique_faces(), fe_q.n_unique_faces());
57 for (unsigned int face_no = 0; face_no < this->n_unique_faces(); ++face_no)
58 this->unit_face_support_points[face_no] =
59 fe_q.get_unit_face_support_points(face_no);
60
61 // initialize unit support points (this makes it possible to assign initial
62 // values to FE_TraceQ). Note that we simply take the points of fe_q but
63 // skip the last ones which are associated with the interior of FE_Q.
64 this->unit_support_points.resize(this->n_dofs_per_cell());
65 for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
66 this->unit_support_points[i] = fe_q.get_unit_support_points()[i];
67
68 // Initialize constraint matrices
69 this->interface_constraints = fe_q.constraints();
70}
71
72
73
74template <int dim, int spacedim>
75std::unique_ptr<FiniteElement<dim, spacedim>>
77{
78 return std::make_unique<FE_TraceQ<dim, spacedim>>(this->degree);
79}
80
81
82
83template <int dim, int spacedim>
84std::string
86{
87 // note that the FETools::get_fe_by_name function depends on the
88 // particular format of the string this function returns, so they have to be
89 // kept in synch
90
91 std::ostringstream namebuf;
92 namebuf << "FE_TraceQ<" << Utilities::dim_string(dim, spacedim) << ">("
93 << this->degree << ")";
94
95 return namebuf.str();
96}
97
98
99
100template <int dim, int spacedim>
101bool
103 const unsigned int shape_index,
104 const unsigned int face_index) const
105{
106 AssertIndexRange(shape_index, this->n_dofs_per_cell());
108
109 // FE_TraceQ shares the numbering of elemental degrees of freedom with FE_Q
110 // except for the missing interior ones (quad dofs in 2d and hex dofs in
111 // 3d). Therefore, it is safe to ask fe_q for the corresponding
112 // information. The assertion 'shape_index < this->n_dofs_per_cell()' will
113 // make sure that we only access the trace dofs.
114 return fe_q.has_support_on_face(shape_index, face_index);
115}
116
117
118
119template <int dim, int spacedim>
120std::pair<Table<2, bool>, std::vector<unsigned int>>
122{
123 Table<2, bool> constant_modes(1, this->n_dofs_per_cell());
124 for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
125 constant_modes(0, i) = true;
126 return std::pair<Table<2, bool>, std::vector<unsigned int>>(
127 constant_modes, std::vector<unsigned int>(1, 0));
128}
129
130template <int dim, int spacedim>
131void
134 const std::vector<Vector<double>> &support_point_values,
135 std::vector<double> & nodal_values) const
136{
137 AssertDimension(support_point_values.size(),
138 this->get_unit_support_points().size());
139 AssertDimension(support_point_values.size(), nodal_values.size());
140 AssertDimension(this->n_dofs_per_cell(), nodal_values.size());
141
142 for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
143 {
144 AssertDimension(support_point_values[i].size(), 1);
145
146 nodal_values[i] = support_point_values[i](0);
147 }
148}
149
150
151template <int dim, int spacedim>
152std::vector<unsigned int>
154{
155 // This constructs FE_TraceQ in exactly the same way as FE_Q except for the
156 // interior degrees of freedom that are not present here (line in 1d, quad
157 // in 2d, hex in 3d).
158 AssertThrow(deg > 0, ExcMessage("FE_TraceQ needs to be of degree > 0."));
159 std::vector<unsigned int> dpo(dim + 1, 1U);
160 dpo[dim] = 0;
161 dpo[0] = 1;
162 for (unsigned int i = 1; i < dim; ++i)
163 dpo[i] = dpo[i - 1] * (deg - 1);
164 return dpo;
165}
166
167
168
169template <int dim, int spacedim>
170bool
172{
173 return fe_q.hp_constraints_are_implemented();
174}
175
176
177template <int dim, int spacedim>
180 const FiniteElement<dim, spacedim> &fe_other,
181 const unsigned int codim) const
182{
183 Assert(codim <= dim, ExcImpossibleInDim(dim));
184 (void)codim;
185
186 // vertex/line/face/cell domination
187 // --------------------------------
188 if (const FE_TraceQ<dim, spacedim> *fe_traceq_other =
189 dynamic_cast<const FE_TraceQ<dim, spacedim> *>(&fe_other))
190 {
191 if (this->degree < fe_traceq_other->degree)
193 else if (this->degree == fe_traceq_other->degree)
195 else
197 }
198 else if (const FE_Nothing<dim> *fe_nothing =
199 dynamic_cast<const FE_Nothing<dim> *>(&fe_other))
200 {
201 if (fe_nothing->is_dominating())
203 else
204 // the FE_Nothing has no degrees of freedom and it is typically used
205 // in a context where we don't require any continuity along the
206 // interface
208 }
209
210 Assert(false, ExcNotImplemented());
212}
213
214
215
216template <int dim, int spacedim>
217void
219 const FiniteElement<dim, spacedim> &source_fe,
220 FullMatrix<double> & interpolation_matrix,
221 const unsigned int face_no) const
222{
223 get_subface_interpolation_matrix(source_fe,
225 interpolation_matrix,
226 face_no);
227}
228
229
230
231template <int dim, int spacedim>
232void
234 const FiniteElement<dim, spacedim> &x_source_fe,
235 const unsigned int subface,
236 FullMatrix<double> & interpolation_matrix,
237 const unsigned int face_no) const
238{
239 // this is the code from FE_FaceQ
240 Assert(interpolation_matrix.n() == this->n_dofs_per_face(face_no),
241 ExcDimensionMismatch(interpolation_matrix.n(),
242 this->n_dofs_per_face(face_no)));
243 Assert(interpolation_matrix.m() == x_source_fe.n_dofs_per_face(face_no),
244 ExcDimensionMismatch(interpolation_matrix.m(),
245 x_source_fe.n_dofs_per_face(face_no)));
246
247 // see if source is a FaceQ element
248 if (const FE_TraceQ<dim, spacedim> *source_fe =
249 dynamic_cast<const FE_TraceQ<dim, spacedim> *>(&x_source_fe))
250 {
251 fe_q.get_subface_interpolation_matrix(source_fe->fe_q,
252 subface,
253 interpolation_matrix,
254 face_no);
255 }
256 else if (dynamic_cast<const FE_Nothing<dim> *>(&x_source_fe) != nullptr)
257 {
258 // nothing to do here, the FE_Nothing has no degrees of freedom anyway
259 }
260 else
262 false,
263 (typename FiniteElement<dim,
264 spacedim>::ExcInterpolationNotImplemented()));
265}
266
267
268
269template <int spacedim>
270FE_TraceQ<1, spacedim>::FE_TraceQ(const unsigned int degree)
271 : FE_FaceQ<1, spacedim>(degree)
272{}
273
274
275
276template <int spacedim>
277std::string
279{
280 // note that the FETools::get_fe_by_name function depends on the
281 // particular format of the string this function returns, so they have to be
282 // kept in synch
283 std::ostringstream namebuf;
284 namebuf << "FE_TraceQ<" << Utilities::dim_string(1, spacedim) << ">("
285 << this->degree << ")";
286
287 return namebuf.str();
288}
289
290
291
292// explicit instantiations
293#include "fe_trace.inst"
294
295
PolynomialType poly_space
FE_Q< dim, spacedim > fe_q
Definition fe_trace.h:141
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
Definition fe_trace.cc:76
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const override
Definition fe_trace.cc:102
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const override
Definition fe_trace.cc:121
virtual bool hp_constraints_are_implemented() const override
Definition fe_trace.cc:171
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final
Definition fe_trace.cc:179
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
Definition fe_trace.cc:133
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix, const unsigned int face_no=0) const override
Definition fe_trace.cc:233
virtual std::string get_name() const override
Definition fe_trace.cc:85
static std::vector< unsigned int > get_dpo_vector(const unsigned int deg)
Definition fe_trace.cc:153
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix, const unsigned int face_no=0) const override
Definition fe_trace.cc:218
FE_TraceQ(unsigned int p)
Definition fe_trace.cc:37
const unsigned int degree
Definition fe_data.h:453
unsigned int n_dofs_per_cell() const
unsigned int n_dofs_per_face(unsigned int face_no=0, unsigned int child=0) const
unsigned int n_unique_faces() const
std::vector< std::vector< Point< dim - 1 > > > unit_face_support_points
Definition fe.h:2447
std::vector< Point< dim > > unit_support_points
Definition fe.h:2440
FullMatrix< double > interface_constraints
Definition fe.h:2428
size_type n() const
size_type m() 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 & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
std::string dim_string(const int dim, const int spacedim)
Definition utilities.cc:556
static const unsigned int invalid_unsigned_int
Definition types.h:213
STL namespace.