Reference documentation for deal.II version GIT 3c37cfefb2 2023-03-24 04:00:03+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_nothing.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2009 - 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 
17 #include <deal.II/fe/fe_nothing.h>
18 
19 #include <memory>
20 
22 
23 
24 template <int dim, int spacedim>
26  const unsigned int n_components,
27  const bool dominate)
28  : FiniteElement<dim, spacedim>(
29  FiniteElementData<dim>(std::vector<unsigned>(dim + 1, 0),
30  type,
31  n_components,
32  0,
33  FiniteElementData<dim>::unknown),
34  std::vector<bool>(),
35  std::vector<ComponentMask>())
36  , dominate(dominate)
37 {
38  Assert(n_components >= 1,
39  ExcMessage("A finite element needs to have at least one "
40  "vector component."));
41 
42  // in most other elements we have to set up all sorts of stuff
43  // here. there isn't much that we have to do here; in particular,
44  // we can simply leave the restriction and prolongation matrices
45  // empty since their proper size is in fact zero given that the
46  // element here has no degrees of freedom
47 }
48 
49 
50 
51 template <int dim, int spacedim>
52 FE_Nothing<dim, spacedim>::FE_Nothing(const unsigned int n_components,
53  const bool dominate)
54  : FE_Nothing<dim, spacedim>(ReferenceCells::get_hypercube<dim>(),
55  n_components,
56  dominate)
57 {}
58 
59 
60 
61 template <int dim, int spacedim>
62 std::unique_ptr<FiniteElement<dim, spacedim>>
64 {
65  return std::make_unique<FE_Nothing<dim, spacedim>>(*this);
66 }
67 
68 
69 
70 template <int dim, int spacedim>
71 std::string
73 {
74  std::ostringstream namebuf;
75  namebuf << "FE_Nothing<" << Utilities::dim_string(dim, spacedim) << ">(";
76 
77  std::vector<std::string> name_components;
78  if (this->reference_cell() != ReferenceCells::get_hypercube<dim>())
79  name_components.push_back(this->reference_cell().to_string());
80  if (this->n_components() > 1)
81  name_components.push_back(std::to_string(this->n_components()));
82  if (dominate)
83  name_components.emplace_back("dominating");
84 
85  for (const std::string &comp : name_components)
86  {
87  namebuf << comp;
88  if (comp != name_components.back())
89  namebuf << ", ";
90  }
91  namebuf << ")";
92 
93  return namebuf.str();
94 }
95 
96 
97 
98 template <int dim, int spacedim>
101 {
102  return flags;
103 }
104 
105 
106 
107 template <int dim, int spacedim>
108 double
109 FE_Nothing<dim, spacedim>::shape_value(const unsigned int /*i*/,
110  const Point<dim> & /*p*/) const
111 {
112  Assert(false, ExcMessage("This element has no shape functions."));
113  return 0;
114 }
115 
116 
117 
118 template <int dim, int spacedim>
119 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
121  const UpdateFlags /*update_flags*/,
122  const Mapping<dim, spacedim> & /*mapping*/,
123  const Quadrature<dim> & /*quadrature*/,
125  spacedim>
126  & /*output_data*/) const
127 {
128  // Create a default data object. Normally we would then
129  // need to resize things to hold the appropriate numbers
130  // of dofs, but in this case all data fields are empty.
131  return std::make_unique<
133 }
134 
135 
136 
137 template <int dim, int spacedim>
138 void
142  const Quadrature<dim> &,
143  const Mapping<dim, spacedim> &,
148  spacedim>
149  &) const
150 {
151  // leave data fields empty
152 }
153 
154 
155 
156 template <int dim, int spacedim>
157 void
160  const unsigned int,
161  const hp::QCollection<dim - 1> &,
162  const Mapping<dim, spacedim> &,
167  spacedim>
168  &) const
169 {
170  // leave data fields empty
171 }
172 
173 
174 
175 template <int dim, int spacedim>
176 void
179  const unsigned int,
180  const unsigned int,
181  const Quadrature<dim - 1> &,
182  const Mapping<dim, spacedim> &,
187  spacedim>
188  &) const
189 {
190  // leave data fields empty
191 }
192 
193 
194 
195 template <int dim, int spacedim>
196 bool
198 {
199  return dominate;
200 }
201 
202 
203 
204 template <int dim, int spacedim>
208  const unsigned int codim) const
209 {
210  Assert(codim <= dim, ExcImpossibleInDim(dim));
211  (void)codim;
212 
213  if (!dominate)
214  // if FE_Nothing does not dominate, there are no requirements
216  else if (dynamic_cast<const FE_Nothing<dim> *>(&fe) != nullptr)
217  // if it does and the other is FE_Nothing, either can dominate
219  else
220  // otherwise we dominate whatever FE is provided
222 }
223 
224 
225 
226 template <int dim, int spacedim>
227 std::vector<std::pair<unsigned int, unsigned int>>
229  const FiniteElement<dim, spacedim> & /*fe_other*/) const
230 {
231  // the FE_Nothing has no
232  // degrees of freedom, so there
233  // are no equivalencies to be
234  // recorded
235  return std::vector<std::pair<unsigned int, unsigned int>>();
236 }
237 
238 
239 template <int dim, int spacedim>
240 std::vector<std::pair<unsigned int, unsigned int>>
242  const FiniteElement<dim, spacedim> & /*fe_other*/) const
243 {
244  // the FE_Nothing has no
245  // degrees of freedom, so there
246  // are no equivalencies to be
247  // recorded
248  return std::vector<std::pair<unsigned int, unsigned int>>();
249 }
250 
251 
252 template <int dim, int spacedim>
253 std::vector<std::pair<unsigned int, unsigned int>>
255  const FiniteElement<dim, spacedim> & /*fe_other*/,
256  const unsigned int) const
257 {
258  // the FE_Nothing has no
259  // degrees of freedom, so there
260  // are no equivalencies to be
261  // recorded
262  return std::vector<std::pair<unsigned int, unsigned int>>();
263 }
264 
265 
266 template <int dim, int spacedim>
267 bool
269 {
270  return true;
271 }
272 
273 
274 
275 template <int dim, int spacedim>
276 void
278  const FiniteElement<dim, spacedim> & /*source_fe*/,
279  FullMatrix<double> &interpolation_matrix) const
280 {
281  // Since this element has no dofs,
282  // the interpolation matrix is necessarily empty.
283  (void)interpolation_matrix;
284 
285  Assert(interpolation_matrix.m() == 0,
286  ExcDimensionMismatch(interpolation_matrix.m(), 0));
287  Assert(interpolation_matrix.n() == 0,
288  ExcDimensionMismatch(interpolation_matrix.n(), 0));
289 }
290 
291 
292 
293 template <int dim, int spacedim>
294 void
296  const FiniteElement<dim, spacedim> & /*source_fe*/,
297  FullMatrix<double> &interpolation_matrix,
298  const unsigned int) const
299 {
300  // since this element has no face dofs, the
301  // interpolation matrix is necessarily empty
302  (void)interpolation_matrix;
303 
304  Assert(interpolation_matrix.m() == 0,
305  ExcDimensionMismatch(interpolation_matrix.m(), 0));
306  Assert(interpolation_matrix.n() == 0,
307  ExcDimensionMismatch(interpolation_matrix.m(), 0));
308 }
309 
310 
311 
312 template <int dim, int spacedim>
313 void
315  const FiniteElement<dim, spacedim> & /*source_fe*/,
316  const unsigned int /*index*/,
317  FullMatrix<double> &interpolation_matrix,
318  const unsigned int) const
319 {
320  // since this element has no face dofs, the
321  // interpolation matrix is necessarily empty
322 
323  (void)interpolation_matrix;
324  Assert(interpolation_matrix.m() == 0,
325  ExcDimensionMismatch(interpolation_matrix.m(), 0));
326  Assert(interpolation_matrix.n() == 0,
327  ExcDimensionMismatch(interpolation_matrix.m(), 0));
328 }
329 
330 
331 
332 template <int dim, int spacedim>
333 std::pair<Table<2, bool>, std::vector<unsigned int>>
335 {
336  // since this element has no dofs, there are no constant modes
337  return {Table<2, bool>{}, std::vector<unsigned int>{}};
338 }
339 
340 
341 
342 // explicit instantiations
343 #include "fe_nothing.inst"
344 
345 
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
Definition: fe_nothing.cc:241
virtual std::string get_name() const override
Definition: fe_nothing.cc:72
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source_fe, const unsigned int index, FullMatrix< double > &interpolation_matrix, const unsigned int face_no=0) const override
Definition: fe_nothing.cc:314
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
Definition: fe_nothing.cc:63
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source_fe, FullMatrix< double > &interpolation_matrix, const unsigned int face_no=0) const override
Definition: fe_nothing.cc:295
FE_Nothing(const ReferenceCell &type, const unsigned int n_components=1, const bool dominate=false)
Definition: fe_nothing.cc:25
virtual void fill_fe_subface_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim - 1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const typename FiniteElement< dim, spacedim >::InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
Definition: fe_nothing.cc:177
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
Definition: fe_nothing.cc:100
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
Definition: fe_nothing.cc:228
virtual double shape_value(const unsigned int i, const Point< dim > &p) const override
Definition: fe_nothing.cc:109
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const hp::QCollection< dim - 1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const typename FiniteElement< dim, spacedim >::InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
Definition: fe_nothing.cc:158
virtual void get_interpolation_matrix(const FiniteElement< dim, spacedim > &source_fe, FullMatrix< double > &interpolation_matrix) const override
Definition: fe_nothing.cc:277
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final
Definition: fe_nothing.cc:206
virtual void fill_fe_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const Quadrature< dim > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const typename FiniteElement< dim, spacedim >::InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
Definition: fe_nothing.cc:139
bool is_dominating() const
Definition: fe_nothing.cc:197
virtual std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
Definition: fe_nothing.cc:120
virtual bool hp_constraints_are_implemented() const override
Definition: fe_nothing.cc:268
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const override
Definition: fe_nothing.cc:334
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_nothing.cc:254
unsigned int n_components() const
size_type n() const
size_type m() const
Definition: point.h:110
#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:1586
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
static ::ExceptionBase & ExcMessage(std::string arg1)
UpdateFlags
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
std::string to_string(const T &t)
Definition: patterns.h:2392
constexpr const ReferenceCell & get_hypercube()
std::string dim_string(const int dim, const int spacedim)
Definition: utilities.cc:556