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_face.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2009 - 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#ifndef dealii_fe_face_h
17#define dealii_fe_face_h
18
19#include <deal.II/base/config.h>
20
23
25
27
28
54template <int dim, int spacedim = dim>
56 : public FE_PolyFace<TensorProductPolynomials<dim - 1>, dim, spacedim>
57{
58public:
64 FE_FaceQ(const unsigned int p);
65
66 virtual std::unique_ptr<FiniteElement<dim, spacedim>>
67 clone() const override;
68
74 virtual std::string
75 get_name() const override;
76
84 virtual void
86 const std::vector<Vector<double>> &support_point_values,
87 std::vector<double> & nodal_values) const override;
88
97 virtual void
99 FullMatrix<double> & matrix,
100 const unsigned int face_no = 0) const override;
101
110 virtual void
112 const FiniteElement<dim, spacedim> &source,
113 const unsigned int subface,
114 FullMatrix<double> & matrix,
115 const unsigned int face_no = 0) const override;
116
121 virtual bool
122 has_support_on_face(const unsigned int shape_index,
123 const unsigned int face_index) const override;
124
147 virtual std::vector<std::pair<unsigned int, unsigned int>>
149 const FiniteElement<dim, spacedim> &fe_other) const override;
150
157 virtual std::vector<std::pair<unsigned int, unsigned int>>
159 const FiniteElement<dim, spacedim> &fe_other) const override;
160
167 virtual std::vector<std::pair<unsigned int, unsigned int>>
169 const unsigned int face_no = 0) const override;
170
175 virtual bool
176 hp_constraints_are_implemented() const override;
177
183 const unsigned int codim = 0) const override final;
184
193 virtual std::pair<Table<2, bool>, std::vector<unsigned int>>
194 get_constant_modes() const override;
195
196private:
200 static std::vector<unsigned int>
201 get_dpo_vector(const unsigned int deg);
202};
203
204
205
216template <int spacedim>
217class FE_FaceQ<1, spacedim> : public FiniteElement<1, spacedim>
218{
219public:
223 FE_FaceQ(const unsigned int p);
224
225 virtual std::unique_ptr<FiniteElement<1, spacedim>>
226 clone() const override;
227
233 virtual std::string
234 get_name() const override;
235
236 // for documentation, see the FiniteElement base class
237 virtual UpdateFlags
238 requires_update_flags(const UpdateFlags update_flags) const override;
239
248 virtual void
250 FullMatrix<double> & matrix,
251 const unsigned int face_no = 0) const override;
252
261 virtual void
263 const FiniteElement<1, spacedim> &source,
264 const unsigned int subface,
265 FullMatrix<double> & matrix,
266 const unsigned int face_no = 0) const override;
267
272 virtual bool
273 has_support_on_face(const unsigned int shape_index,
274 const unsigned int face_index) const override;
275
280 virtual bool
281 hp_constraints_are_implemented() const override;
282
300 virtual std::vector<std::pair<unsigned int, unsigned int>>
302 const FiniteElement<1, spacedim> &fe_other) const override;
303
310 virtual std::vector<std::pair<unsigned int, unsigned int>>
312 const FiniteElement<1, spacedim> &fe_other) const override;
313
320 virtual std::vector<std::pair<unsigned int, unsigned int>>
322 const unsigned int face_no = 0) const override;
323
328 virtual std::pair<Table<2, bool>, std::vector<unsigned int>>
329 get_constant_modes() const override;
330
331protected:
332 /*
333 * NOTE: The following functions have their definitions inlined into the class
334 * declaration because we otherwise run into a compiler error with MS Visual
335 * Studio.
336 */
337
338
339 virtual std::unique_ptr<typename FiniteElement<1, spacedim>::InternalDataBase>
341 const UpdateFlags /*update_flags*/,
342 const Mapping<1, spacedim> & /*mapping*/,
343 const Quadrature<1> & /*quadrature*/,
345 spacedim>
346 & /*output_data*/) const override
347 {
348 return std::make_unique<
350 }
351
352 using FiniteElement<1, spacedim>::get_face_data;
353
354 std::unique_ptr<typename FiniteElement<1, spacedim>::InternalDataBase>
356 const UpdateFlags update_flags,
357 const Mapping<1, spacedim> & /*mapping*/,
358 const hp::QCollection<0> &quadrature,
360 spacedim>
361 & /*output_data*/) const override
362 {
363 AssertDimension(quadrature.size(), 1);
364
365 // generate a new data object and initialize some fields
366 auto data_ptr =
367 std::make_unique<typename FiniteElement<1, spacedim>::InternalDataBase>();
368 data_ptr->update_each = requires_update_flags(update_flags);
369
370 const unsigned int n_q_points = quadrature[0].size();
371 AssertDimension(n_q_points, 1);
372 (void)n_q_points;
373
374 // No derivatives of this element are implemented.
375 if (data_ptr->update_each & update_gradients ||
376 data_ptr->update_each & update_hessians)
377 {
378 Assert(false, ExcNotImplemented());
379 }
380
381 return data_ptr;
382 }
383
384 std::unique_ptr<typename FiniteElement<1, spacedim>::InternalDataBase>
386 const UpdateFlags update_flags,
387 const Mapping<1, spacedim> &mapping,
388 const Quadrature<0> & quadrature,
390 spacedim>
391 &output_data) const override
392 {
393 return get_face_data(update_flags,
394 mapping,
395 hp::QCollection<0>(quadrature),
396 output_data);
397 }
398
399 virtual void
402 const CellSimilarity::Similarity cell_similarity,
403 const Quadrature<1> & quadrature,
404 const Mapping<1, spacedim> & mapping,
405 const typename Mapping<1, spacedim>::InternalDataBase & mapping_internal,
407 & mapping_data,
408 const typename FiniteElement<1, spacedim>::InternalDataBase &fe_internal,
410 spacedim>
411 &output_data) const override;
412
413 using FiniteElement<1, spacedim>::fill_fe_face_values;
414
415 virtual void
418 const unsigned int face_no,
419 const hp::QCollection<0> & quadrature,
420 const Mapping<1, spacedim> & mapping,
421 const typename Mapping<1, spacedim>::InternalDataBase & mapping_internal,
423 & mapping_data,
424 const typename FiniteElement<1, spacedim>::InternalDataBase &fe_internal,
426 spacedim>
427 &output_data) const override;
428
429 virtual void
432 const unsigned int face_no,
433 const unsigned int sub_no,
434 const Quadrature<0> & quadrature,
435 const Mapping<1, spacedim> & mapping,
436 const typename Mapping<1, spacedim>::InternalDataBase & mapping_internal,
438 & mapping_data,
439 const typename FiniteElement<1, spacedim>::InternalDataBase &fe_internal,
441 spacedim>
442 &output_data) const override;
443
444private:
448 static std::vector<unsigned int>
449 get_dpo_vector(const unsigned int deg);
450};
451
452
453
476template <int dim, int spacedim = dim>
477class FE_FaceP : public FE_PolyFace<PolynomialSpace<dim - 1>, dim, spacedim>
478{
479public:
485 FE_FaceP(unsigned int p);
486
487 virtual std::unique_ptr<FiniteElement<dim, spacedim>>
488 clone() const override;
489
495 virtual std::string
496 get_name() const override;
497
506 virtual void
508 FullMatrix<double> & matrix,
509 const unsigned int face_no = 0) const override;
510
519 virtual void
521 const FiniteElement<dim, spacedim> &source,
522 const unsigned int subface,
523 FullMatrix<double> & matrix,
524 const unsigned int face_no = 0) const override;
525
530 virtual bool
531 has_support_on_face(const unsigned int shape_index,
532 const unsigned int face_index) const override;
533
538 virtual bool
539 hp_constraints_are_implemented() const override;
540
546 const unsigned int codim = 0) const override final;
547
554 virtual std::pair<Table<2, bool>, std::vector<unsigned int>>
555 get_constant_modes() const override;
556
557private:
561 static std::vector<unsigned int>
562 get_dpo_vector(const unsigned int deg);
563};
564
565
566
571template <int spacedim>
572class FE_FaceP<1, spacedim> : public FE_FaceQ<1, spacedim>
573{
574public:
578 FE_FaceP(const unsigned int p);
579
583 std::string
584 get_name() const override;
585};
586
587
589
590#endif
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const override
Definition fe_face.cc:994
virtual bool hp_constraints_are_implemented() const override
Definition fe_face.cc:824
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const override
Definition fe_face.cc:797
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final
Definition fe_face.cc:833
virtual std::string get_name() const override
Definition fe_face.cc:781
static std::vector< unsigned int > get_dpo_vector(const unsigned int deg)
Definition fe_face.cc:808
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_face.cc:887
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
Definition fe_face.cc:772
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix, const unsigned int face_no=0) const override
Definition fe_face.cc:872
std::unique_ptr< typename FiniteElement< 1, spacedim >::InternalDataBase > get_face_data(const UpdateFlags update_flags, const Mapping< 1, spacedim > &, const hp::QCollection< 0 > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< 1, spacedim > &) const override
Definition fe_face.h:355
virtual std::unique_ptr< typename FiniteElement< 1, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Mapping< 1, spacedim > &, const Quadrature< 1 > &, ::internal::FEValuesImplementation::FiniteElementRelatedData< 1, spacedim > &) const override
Definition fe_face.h:340
std::unique_ptr< typename FiniteElement< 1, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags update_flags, const Mapping< 1, spacedim > &mapping, const Quadrature< 0 > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< 1, spacedim > &output_data) const override
Definition fe_face.h:385
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final
Definition fe_face.cc:448
virtual std::string get_name() const override
Definition fe_face.cc:131
static std::vector< unsigned int > get_dpo_vector(const unsigned int deg)
Definition fe_face.cc:266
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_face.cc:496
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const override
Definition fe_face.cc:485
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix, const unsigned int face_no=0) const override
Definition fe_face.cc:147
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
Definition fe_face.cc:299
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
Definition fe_face.cc:122
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_face.cc:162
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const override
Definition fe_face.cc:255
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
Definition fe_face.cc:288
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_face.cc:369
virtual bool hp_constraints_are_implemented() const override
Definition fe_face.cc:279
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
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > get_face_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const hp::QCollection< dim - 1 > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
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
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
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
Abstract base class for mapping classes.
Definition mapping.h:317
unsigned int size() const
Definition collection.h:265
#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)
#define AssertDimension(dim1, dim2)
UpdateFlags
@ update_hessians
Second derivatives of shape functions.
@ update_gradients
Shape function gradients.