Reference documentation for deal.II version 9.4.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\}}\)
fe_q.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2000 - 2021 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
18
19#include <deal.II/fe/fe_dgq.h>
22#include <deal.II/fe/fe_q.h>
26
27#include <deal.II/lac/vector.h>
28
29#include <memory>
30#include <sstream>
31#include <vector>
32
34
35
36namespace internal
37{
38 namespace FE_Q
39 {
40 namespace
41 {
42 std::vector<Point<1>>
43 get_QGaussLobatto_points(const unsigned int degree)
44 {
45 if (degree > 0)
46 return QGaussLobatto<1>(degree + 1).get_points();
47 else
48 {
49 using FEQ = ::FE_Q_Base<1, 1>;
50 AssertThrow(false, FEQ::ExcFEQCannotHaveDegree0());
51 }
52 return std::vector<Point<1>>();
53 }
54 } // namespace
55 } // namespace FE_Q
56} // namespace internal
57
58
59
60template <int dim, int spacedim>
61FE_Q<dim, spacedim>::FE_Q(const unsigned int degree)
62 : FE_Q_Base<dim, spacedim>(
65 internal::FE_Q::get_QGaussLobatto_points(degree))),
66 FiniteElementData<dim>(this->get_dpo_vector(degree),
67 1,
68 degree,
69 FiniteElementData<dim>::H1),
70 std::vector<bool>(1, false))
71{
72 this->initialize(internal::FE_Q::get_QGaussLobatto_points(degree));
73}
74
75
76
77template <int dim, int spacedim>
79 : FE_Q_Base<dim, spacedim>(
81 Polynomials::generate_complete_Lagrange_basis(points.get_points())),
82 FiniteElementData<dim>(this->get_dpo_vector(points.size() - 1),
83 1,
84 points.size() - 1,
85 FiniteElementData<dim>::H1),
86 std::vector<bool>(1, false))
87{
88 this->initialize(points.get_points());
89}
90
91
92
93template <int dim, int spacedim>
94std::string
96{
97 // note that the FETools::get_fe_by_name function depends on the
98 // particular format of the string this function returns, so they have to be
99 // kept in synch
100
101 std::ostringstream namebuf;
102 bool equidistant = true;
103 std::vector<double> points(this->degree + 1);
104
105 // Decode the support points in one coordinate direction.
106 TensorProductPolynomials<dim> *poly_space_derived_ptr =
107 dynamic_cast<TensorProductPolynomials<dim> *>(this->poly_space.get());
108 std::vector<unsigned int> lexicographic =
109 poly_space_derived_ptr->get_numbering_inverse();
110 for (unsigned int j = 0; j <= this->degree; ++j)
111 points[j] = this->unit_support_points[lexicographic[j]][0];
112
113 // Check whether the support points are equidistant.
114 for (unsigned int j = 0; j <= this->degree; ++j)
115 if (std::fabs(points[j] - static_cast<double>(j) / this->degree) > 1e-15)
116 {
117 equidistant = false;
118 break;
119 }
120
121 if (equidistant == true)
122 {
123 if (this->degree > 2)
124 namebuf << "FE_Q<" << Utilities::dim_string(dim, spacedim)
125 << ">(QIterated(QTrapezoid()," << this->degree << "))";
126 else
127 namebuf << "FE_Q<" << Utilities::dim_string(dim, spacedim) << ">("
128 << this->degree << ")";
129 }
130 else
131 {
132 // Check whether the support points come from QGaussLobatto.
133 const QGaussLobatto<1> points_gl(this->degree + 1);
134 bool gauss_lobatto = true;
135 for (unsigned int j = 0; j <= this->degree; ++j)
136 if (points[j] != points_gl.point(j)(0))
137 {
138 gauss_lobatto = false;
139 break;
140 }
141 if (gauss_lobatto == true)
142 namebuf << "FE_Q<" << Utilities::dim_string(dim, spacedim) << ">("
143 << this->degree << ")";
144 else
145 namebuf << "FE_Q<" << Utilities::dim_string(dim, spacedim)
146 << ">(QUnknownNodes(" << this->degree << "))";
147 }
148 return namebuf.str();
149}
150
151
152
153template <int dim, int spacedim>
154void
156 const std::vector<Vector<double>> &support_point_values,
157 std::vector<double> & nodal_values) const
158{
159 AssertDimension(support_point_values.size(),
160 this->get_unit_support_points().size());
161 AssertDimension(support_point_values.size(), nodal_values.size());
162 AssertDimension(this->n_dofs_per_cell(), nodal_values.size());
163
164 for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
165 {
166 AssertDimension(support_point_values[i].size(), 1);
167
168 nodal_values[i] = support_point_values[i](0);
169 }
170}
171
172
173
174template <int dim, int spacedim>
175std::unique_ptr<FiniteElement<dim, spacedim>>
177{
178 return std::make_unique<FE_Q<dim, spacedim>>(*this);
179}
180
181
182
183template <int dim, int spacedim>
186 const FiniteElement<dim, spacedim> &fe_other,
187 const unsigned int codim) const
188{
189 Assert(codim <= dim, ExcImpossibleInDim(dim));
190
191 // vertex/line/face domination
192 // (if fe_other is derived from FE_DGQ)
193 // ------------------------------------
194 if (codim > 0)
195 if (dynamic_cast<const FE_DGQ<dim, spacedim> *>(&fe_other) != nullptr)
196 // there are no requirements between continuous and discontinuous elements
198
199 // vertex/line/face domination
200 // (if fe_other is not derived from FE_DGQ)
201 // & cell domination
202 // ----------------------------------------
203 if (const FE_Q<dim, spacedim> *fe_q_other =
204 dynamic_cast<const FE_Q<dim, spacedim> *>(&fe_other))
205 {
206 if (this->degree < fe_q_other->degree)
208 else if (this->degree == fe_q_other->degree)
210 else
212 }
213 else if (const FE_SimplexP<dim, spacedim> *fe_p_other =
214 dynamic_cast<const FE_SimplexP<dim, spacedim> *>(&fe_other))
215 {
216 if (this->degree < fe_p_other->degree)
218 else if (this->degree == fe_p_other->degree)
220 else
222 }
223 else if (const FE_WedgeP<dim, spacedim> *fe_wp_other =
224 dynamic_cast<const FE_WedgeP<dim, spacedim> *>(&fe_other))
225 {
226 if (this->degree < fe_wp_other->degree)
228 else if (this->degree == fe_wp_other->degree)
230 else
232 }
233 else if (const FE_PyramidP<dim, spacedim> *fe_pp_other =
234 dynamic_cast<const FE_PyramidP<dim, spacedim> *>(&fe_other))
235 {
236 if (this->degree < fe_pp_other->degree)
238 else if (this->degree == fe_pp_other->degree)
240 else
242 }
243 else if (const FE_Nothing<dim, spacedim> *fe_nothing =
244 dynamic_cast<const FE_Nothing<dim, spacedim> *>(&fe_other))
245 {
246 if (fe_nothing->is_dominating())
248 else
249 // the FE_Nothing has no degrees of freedom and it is typically used
250 // in a context where we don't require any continuity along the
251 // interface
253 }
254
255 Assert(false, ExcNotImplemented());
257}
258
259
260// explicit instantiations
261#include "fe_q.inst"
262
Definition: fe_dgq.h:111
void initialize(const std::vector< Point< 1 > > &support_points_1d)
Definition: fe_q.h:549
FE_Q(const unsigned int p)
Definition: fe_q.cc:61
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_q.cc:155
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
Definition: fe_q.cc:176
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final
Definition: fe_q.cc:185
virtual std::string get_name() const override
Definition: fe_q.cc:95
const unsigned int degree
Definition: fe_data.h:449
const Point< dim > & point(const unsigned int i) const
const std::vector< Point< dim > > & get_points() const
const std::vector< unsigned int > & get_numbering_inverse() const
Definition: vector.h:109
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
Definition: exceptions.h:1473
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1667
#define AssertThrow(cond, exc)
Definition: exceptions.h:1583
Expression fabs(const Expression &x)
std::vector< unsigned int > get_dpo_vector(const unsigned int degree)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
std::vector< Polynomial< double > > generate_complete_Lagrange_basis(const std::vector< Point< 1 > > &points)
Definition: polynomial.cc:702
std::string dim_string(const int dim, const int spacedim)
Definition: utilities.cc:558
std::vector< Point< dim > > unit_support_points(const std::vector< Point< 1 > > &line_support_points, const std::vector< unsigned int > &renumbering)
STL namespace.