Reference documentation for deal.II version GIT a3f17f8a20 2023-06-02 10:50: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\}}\)
fe_dgp.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2002 - 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 
16 
17 #include <deal.II/fe/fe_dgp.h>
18 #include <deal.II/fe/fe_nothing.h>
19 #include <deal.II/fe/fe_tools.h>
20 
21 #include <memory>
22 #include <sstream>
23 
24 
26 
27 template <int dim, int spacedim>
28 FE_DGP<dim, spacedim>::FE_DGP(const unsigned int degree)
29  : FE_Poly<dim, spacedim>(
30  PolynomialSpace<dim>(
31  Polynomials::Legendre::generate_complete_basis(degree)),
32  FiniteElementData<dim>(get_dpo_vector(degree),
33  1,
34  degree,
35  FiniteElementData<dim>::L2),
36  std::vector<bool>(
37  FiniteElementData<dim>(get_dpo_vector(degree), 1, degree)
38  .n_dofs_per_cell(),
39  true),
40  std::vector<ComponentMask>(
41  FiniteElementData<dim>(get_dpo_vector(degree), 1, degree)
42  .n_dofs_per_cell(),
43  std::vector<bool>(1, true)))
44 {
45  // Reinit the vectors of restriction and prolongation matrices to the right
46  // sizes
48  // Fill prolongation matrices with embedding operators
49  if (dim == spacedim)
50  {
52  // Fill restriction matrices with L2-projection
54  }
55 }
56 
57 
58 template <int dim, int spacedim>
59 std::string
61 {
62  // note that the FETools::get_fe_by_name function depends on the
63  // particular format of the string this function returns, so they have to be
64  // kept in sync
65 
66  std::ostringstream namebuf;
67  namebuf << "FE_DGP<" << Utilities::dim_string(dim, spacedim) << ">("
68  << this->degree << ")";
69 
70  return namebuf.str();
71 }
72 
73 
74 
75 template <int dim, int spacedim>
76 std::unique_ptr<FiniteElement<dim, spacedim>>
78 {
79  return std::make_unique<FE_DGP<dim, spacedim>>(*this);
80 }
81 
82 
83 
84 //---------------------------------------------------------------------------
85 // Auxiliary functions
86 //---------------------------------------------------------------------------
87 
88 
89 template <int dim, int spacedim>
90 std::vector<unsigned int>
92 {
93  std::vector<unsigned int> dpo(dim + 1, 0U);
94  dpo[dim] = deg + 1;
95  for (unsigned int i = 1; i < dim; ++i)
96  {
97  dpo[dim] *= deg + 1 + i;
98  dpo[dim] /= i + 1;
99  }
100  return dpo;
101 }
102 
103 
104 
105 template <int dim, int spacedim>
106 void
108  const FiniteElement<dim, spacedim> &x_source_fe,
109  FullMatrix<double> & interpolation_matrix,
110  const unsigned int) const
111 {
112  // this is only implemented, if the source FE is also a DGP element. in that
113  // case, both elements have no dofs on their faces and the face
114  // interpolation matrix is necessarily empty -- i.e. there isn't much we
115  // need to do here.
116  (void)interpolation_matrix;
117  using FE = FiniteElement<dim, spacedim>;
118  using FEDGP = FE_DGP<dim, spacedim>;
119  AssertThrow((x_source_fe.get_name().find("FE_DGP<") == 0) ||
120  (dynamic_cast<const FEDGP *>(&x_source_fe) != nullptr),
121  typename FE::ExcInterpolationNotImplemented());
122 
123  Assert(interpolation_matrix.m() == 0,
124  ExcDimensionMismatch(interpolation_matrix.m(), 0));
125  Assert(interpolation_matrix.n() == 0,
126  ExcDimensionMismatch(interpolation_matrix.n(), 0));
127 }
128 
129 
130 
131 template <int dim, int spacedim>
132 void
134  const FiniteElement<dim, spacedim> &x_source_fe,
135  const unsigned int,
136  FullMatrix<double> &interpolation_matrix,
137  const unsigned int) const
138 {
139  // this is only implemented, if the source FE is also a DGP element. in that
140  // case, both elements have no dofs on their faces and the face
141  // interpolation matrix is necessarily empty -- i.e. there isn't much we
142  // need to do here.
143  (void)interpolation_matrix;
144  using FE = FiniteElement<dim, spacedim>;
145  using FEDGP = FE_DGP<dim, spacedim>;
146  AssertThrow((x_source_fe.get_name().find("FE_DGP<") == 0) ||
147  (dynamic_cast<const FEDGP *>(&x_source_fe) != nullptr),
148  typename FE::ExcInterpolationNotImplemented());
149 
150  Assert(interpolation_matrix.m() == 0,
151  ExcDimensionMismatch(interpolation_matrix.m(), 0));
152  Assert(interpolation_matrix.n() == 0,
153  ExcDimensionMismatch(interpolation_matrix.n(), 0));
154 }
155 
156 
157 
158 template <int dim, int spacedim>
159 bool
161 {
162  return true;
163 }
164 
165 
166 
167 template <int dim, int spacedim>
168 std::vector<std::pair<unsigned int, unsigned int>>
170  const FiniteElement<dim, spacedim> &fe_other) const
171 {
172  // there are no such constraints for DGP elements at all
173  if (dynamic_cast<const FE_DGP<dim, spacedim> *>(&fe_other) != nullptr)
174  return std::vector<std::pair<unsigned int, unsigned int>>();
175  else
176  {
177  Assert(false, ExcNotImplemented());
178  return std::vector<std::pair<unsigned int, unsigned int>>();
179  }
180 }
181 
182 
183 
184 template <int dim, int spacedim>
185 std::vector<std::pair<unsigned int, unsigned int>>
187  const FiniteElement<dim, spacedim> &fe_other) const
188 {
189  // there are no such constraints for DGP elements at all
190  if (dynamic_cast<const FE_DGP<dim, spacedim> *>(&fe_other) != nullptr)
191  return std::vector<std::pair<unsigned int, unsigned int>>();
192  else
193  {
194  Assert(false, ExcNotImplemented());
195  return std::vector<std::pair<unsigned int, unsigned int>>();
196  }
197 }
198 
199 
200 
201 template <int dim, int spacedim>
202 std::vector<std::pair<unsigned int, unsigned int>>
204  const FiniteElement<dim, spacedim> &fe_other,
205  const unsigned int) const
206 {
207  // there are no such constraints for DGP elements at all
208  if (dynamic_cast<const FE_DGP<dim, spacedim> *>(&fe_other) != nullptr)
209  return std::vector<std::pair<unsigned int, unsigned int>>();
210  else
211  {
212  Assert(false, ExcNotImplemented());
213  return std::vector<std::pair<unsigned int, unsigned int>>();
214  }
215 }
216 
217 
218 
219 template <int dim, int spacedim>
222  const FiniteElement<dim, spacedim> &fe_other,
223  const unsigned int codim) const
224 {
225  Assert(codim <= dim, ExcImpossibleInDim(dim));
226 
227  // vertex/line/face domination
228  // ---------------------------
229  if (codim > 0)
230  // this is a discontinuous element, so by definition there will
231  // be no constraints wherever this element comes together with
232  // any other kind of element
234 
235  // cell domination
236  // ---------------
237  if (const FE_DGP<dim, spacedim> *fe_dgp_other =
238  dynamic_cast<const FE_DGP<dim, spacedim> *>(&fe_other))
239  {
240  if (this->degree < fe_dgp_other->degree)
242  else if (this->degree == fe_dgp_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 
265 template <int dim, int spacedim>
266 bool
268  const unsigned int) const
269 {
270  // all shape functions have support on all faces
271  return true;
272 }
273 
274 
275 
276 template <int dim, int spacedim>
277 std::pair<Table<2, bool>, std::vector<unsigned int>>
279 {
280  Table<2, bool> constant_modes(1, this->n_dofs_per_cell());
281  constant_modes(0, 0) = true;
282  return std::pair<Table<2, bool>, std::vector<unsigned int>>(
283  constant_modes, std::vector<unsigned int>(1, 0));
284 }
285 
286 
287 
288 template <int dim, int spacedim>
289 std::size_t
291 {
292  Assert(false, ExcNotImplemented());
293  return 0;
294 }
295 
296 
297 
298 // explicit instantiations
299 #include "fe_dgp.inst"
300 
301 
Definition: fe_dgp.h:312
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
Definition: fe_dgp.cc:186
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const override
Definition: fe_dgp.cc:267
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const override
Definition: fe_dgp.cc:278
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_dgp.cc:133
virtual std::string get_name() const override
Definition: fe_dgp.cc:60
static std::vector< unsigned int > get_dpo_vector(const unsigned int degree)
Definition: fe_dgp.cc:91
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
Definition: fe_dgp.cc:169
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix, const unsigned int face_no=0) const override
Definition: fe_dgp.cc:107
FE_DGP(const unsigned int p)
Definition: fe_dgp.cc:28
virtual std::size_t memory_consumption() const override
Definition: fe_dgp.cc:290
virtual 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
Definition: fe_dgp.cc:203
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
Definition: fe_dgp.cc:77
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final
Definition: fe_dgp.cc:221
virtual bool hp_constraints_are_implemented() const override
Definition: fe_dgp.cc:160
virtual std::string get_name() const =0
std::vector< std::vector< FullMatrix< double > > > restriction
Definition: fe.h:2402
void reinit_restriction_and_prolongation_matrices(const bool isotropic_restriction_only=false, const bool isotropic_prolongation_only=false)
std::vector< std::vector< FullMatrix< double > > > prolongation
Definition: fe.h:2416
size_type n() const
size_type m() const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:474
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:475
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
Definition: exceptions.h:1614
static ::ExceptionBase & ExcNotImplemented()
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1703
void compute_embedding_matrices(const FiniteElement< dim, spacedim > &fe, std::vector< std::vector< FullMatrix< number >>> &matrices, const bool isotropic_only=false, const double threshold=1.e-12)
void compute_projection_matrices(const FiniteElement< dim, spacedim > &fe, std::vector< std::vector< FullMatrix< number >>> &matrices, const bool isotropic_only=false)
std::vector< unsigned int > get_dpo_vector(const unsigned int degree)
static const char U
void L2(Vector< number > &result, const FEValuesBase< dim > &fe, const std::vector< double > &input, const double factor=1.)
Definition: l2.h:160
std::string dim_string(const int dim, const int spacedim)
Definition: utilities.cc:556