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_poly_face.h
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#ifndef dealii_fe_poly_face_h
17#define dealii_fe_poly_face_h
18
19
20#include <deal.II/base/config.h>
21
23
24#include <deal.II/fe/fe.h>
25
26#include <memory>
27
29
58template <class PolynomialType,
59 int dim = PolynomialType::dimension + 1,
60 int spacedim = dim>
61class FE_PolyFace : public FiniteElement<dim, spacedim>
62{
63public:
67 FE_PolyFace(const PolynomialType & poly_space,
68 const FiniteElementData<dim> &fe_data,
69 const std::vector<bool> & restriction_is_additive_flags);
70
75 unsigned int
76 get_degree() const;
77
78 // for documentation, see the FiniteElement base class
79 virtual UpdateFlags
80 requires_update_flags(const UpdateFlags update_flags) const override;
81
82protected:
83 /*
84 * NOTE: The following functions have their definitions inlined into the class
85 * declaration because we otherwise run into a compiler error with MS Visual
86 * Studio.
87 */
88
89
90 virtual std::unique_ptr<
93 const UpdateFlags update_flags,
94 const Mapping<dim, spacedim> &mapping,
95 const Quadrature<dim> & quadrature,
97 spacedim>
98 &output_data) const override
99 {
100 (void)update_flags;
101 (void)mapping;
102 (void)quadrature;
103 (void)output_data;
104 return std::make_unique<InternalData>();
105 }
106
107 using FiniteElement<dim, spacedim>::get_face_data;
108
109 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
111 const UpdateFlags update_flags,
112 const Mapping<dim, spacedim> & mapping,
113 const hp::QCollection<dim - 1> &quadrature,
115 spacedim>
116 &output_data) const override
117 {
118 (void)mapping;
119 (void)output_data;
120 AssertDimension(quadrature.size(), 1);
121
122 // generate a new data object and
123 // initialize some fields
124 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
125 data_ptr = std::make_unique<InternalData>();
126 auto &data = dynamic_cast<InternalData &>(*data_ptr);
127 data.update_each = requires_update_flags(update_flags);
128
129 const unsigned int n_q_points = quadrature[0].size();
130
131 // some scratch arrays
132 std::vector<double> values(0);
133 std::vector<Tensor<1, dim - 1>> grads(0);
134 std::vector<Tensor<2, dim - 1>> grad_grads(0);
135 std::vector<Tensor<3, dim - 1>> empty_vector_of_3rd_order_tensors;
136 std::vector<Tensor<4, dim - 1>> empty_vector_of_4th_order_tensors;
137
138 // initialize fields only if really
139 // necessary. otherwise, don't
140 // allocate memory
141 if (data.update_each & update_values)
142 {
143 values.resize(poly_space.n());
144 data.shape_values.resize(poly_space.n(),
145 std::vector<double>(n_q_points));
146 for (unsigned int i = 0; i < n_q_points; ++i)
147 {
148 poly_space.evaluate(quadrature[0].point(i),
149 values,
150 grads,
151 grad_grads,
152 empty_vector_of_3rd_order_tensors,
153 empty_vector_of_4th_order_tensors);
154
155 for (unsigned int k = 0; k < poly_space.n(); ++k)
156 data.shape_values[k][i] = values[k];
157 }
158 }
159 // No derivatives of this element
160 // are implemented.
161 if (data.update_each & update_gradients ||
162 data.update_each & update_hessians)
163 {
164 Assert(false, ExcNotImplemented());
165 }
166
167 return data_ptr;
168 }
169
170 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
172 const UpdateFlags update_flags,
173 const Mapping<dim, spacedim> &mapping,
174 const Quadrature<dim - 1> & quadrature,
176 spacedim>
177 &output_data) const override
178 {
179 return get_face_data(
180 update_flags,
181 mapping,
183 ReferenceCells::get_hypercube<dim - 1>(), quadrature)),
184 output_data);
185 }
186
187 virtual void
190 const CellSimilarity::Similarity cell_similarity,
191 const Quadrature<dim> & quadrature,
192 const Mapping<dim, spacedim> & mapping,
193 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
195 & mapping_data,
196 const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
198 spacedim>
199 &output_data) const override;
200
201 using FiniteElement<dim, spacedim>::fill_fe_face_values;
202
203 virtual void
206 const unsigned int face_no,
207 const hp::QCollection<dim - 1> & quadrature,
208 const Mapping<dim, spacedim> & mapping,
209 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
211 & mapping_data,
212 const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
214 spacedim>
215 &output_data) const override;
216
217 virtual void
220 const unsigned int face_no,
221 const unsigned int sub_no,
222 const Quadrature<dim - 1> & quadrature,
223 const Mapping<dim, spacedim> & mapping,
224 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
226 & mapping_data,
227 const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
229 spacedim>
230 &output_data) const override;
231
238 class InternalData : public FiniteElement<dim, spacedim>::InternalDataBase
239 {
240 public:
254 std::vector<std::vector<double>> shape_values;
255 };
256
261 PolynomialType poly_space;
262};
263
267
268#endif
std::vector< std::vector< double > > shape_values
PolynomialType poly_space
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim - 1 > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
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
FE_PolyFace(const PolynomialType &poly_space, const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags)
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
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
unsigned int get_degree() const
const std::vector< bool > restriction_is_additive_flags
Definition fe.h:2570
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_values
Shape function values.
@ update_gradients
Shape function gradients.