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.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2004 - 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_h
17#define dealii_fe_poly_h
18
19
20#include <deal.II/base/config.h>
21
24
25#include <deal.II/fe/fe.h>
26
27#include <memory>
28
30
75template <int dim, int spacedim = dim>
76class FE_Poly : public FiniteElement<dim, spacedim>
77{
78public:
83 const FiniteElementData<dim> & fe_data,
84 const std::vector<bool> & restriction_is_additive_flags,
85 const std::vector<ComponentMask> &nonzero_components);
86
90 FE_Poly(const FE_Poly &fe);
91
96 unsigned int
97 get_degree() const;
98
99 // for documentation, see the FiniteElement base class
100 virtual UpdateFlags
101 requires_update_flags(const UpdateFlags update_flags) const override;
102
108
118 std::vector<unsigned int>
120
127 std::vector<unsigned int>
129
135 virtual double
136 shape_value(const unsigned int i, const Point<dim> &p) const override;
137
148 virtual double
149 shape_value_component(const unsigned int i,
150 const Point<dim> & p,
151 const unsigned int component) const override;
152
158 virtual Tensor<1, dim>
159 shape_grad(const unsigned int i, const Point<dim> &p) const override;
160
171 virtual Tensor<1, dim>
172 shape_grad_component(const unsigned int i,
173 const Point<dim> & p,
174 const unsigned int component) const override;
175
181 virtual Tensor<2, dim>
182 shape_grad_grad(const unsigned int i, const Point<dim> &p) const override;
183
194 virtual Tensor<2, dim>
195 shape_grad_grad_component(const unsigned int i,
196 const Point<dim> & p,
197 const unsigned int component) const override;
198
204 virtual Tensor<3, dim>
205 shape_3rd_derivative(const unsigned int i,
206 const Point<dim> & p) const override;
207
218 virtual Tensor<3, dim>
219 shape_3rd_derivative_component(const unsigned int i,
220 const Point<dim> & p,
221 const unsigned int component) const override;
222
228 virtual Tensor<4, dim>
229 shape_4th_derivative(const unsigned int i,
230 const Point<dim> & p) const override;
231
242 virtual Tensor<4, dim>
243 shape_4th_derivative_component(const unsigned int i,
244 const Point<dim> & p,
245 const unsigned int component) const override;
246
250 virtual std::size_t
251 memory_consumption() const override;
252
253protected:
254 /*
255 * NOTE: The following function has its definition inlined into the class
256 * declaration because we otherwise run into a compiler error with MS Visual
257 * Studio.
258 */
259
260
261 virtual std::unique_ptr<
264 const UpdateFlags update_flags,
265 const Mapping<dim, spacedim> &mapping,
266 const Quadrature<dim> & quadrature,
268 spacedim>
269 &output_data) const override
270 {
271 (void)mapping;
272
273 // generate a new data object and
274 // initialize some fields
275 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
276 data_ptr = std::make_unique<InternalData>();
277 auto &data = dynamic_cast<InternalData &>(*data_ptr);
278 data.update_each = requires_update_flags(update_flags);
279
280 const unsigned int n_q_points = quadrature.size();
281
282 // initialize some scratch arrays. we need them for the underlying
283 // polynomial to put the values and derivatives of shape functions
284 // to put there, depending on what the user requested
285 std::vector<double> values(
286 update_flags & update_values ? this->n_dofs_per_cell() : 0);
287 std::vector<Tensor<1, dim>> grads(
288 update_flags & update_gradients ? this->n_dofs_per_cell() : 0);
289 std::vector<Tensor<2, dim>> grad_grads(
290 update_flags & update_hessians ? this->n_dofs_per_cell() : 0);
291 std::vector<Tensor<3, dim>> third_derivatives(
292 update_flags & update_3rd_derivatives ? this->n_dofs_per_cell() : 0);
293 std::vector<Tensor<4, dim>>
294 fourth_derivatives; // won't be needed, so leave empty
295
296 // now also initialize fields the fields of this class's own
297 // temporary storage, depending on what we need for the given
298 // update flags.
299 //
300 // there is one exception from the rule: if we are dealing with
301 // cells (i.e., if this function is not called via
302 // get_(sub)face_data()), then we can already store things in the
303 // final location where FEValues::reinit() later wants to see
304 // things. we then don't need the intermediate space. we determine
305 // whether we are on a cell by asking whether the number of
306 // elements in the output array equals the number of quadrature
307 // points (yes, it's a cell) or not (because in that case the
308 // number of quadrature points we use here equals the number of
309 // quadrature points summed over *all* faces or subfaces, whereas
310 // the number of output slots equals the number of quadrature
311 // points on only *one* face)
312 if ((update_flags & update_values) &&
313 !((output_data.shape_values.n_rows() > 0) &&
314 (output_data.shape_values.n_cols() == n_q_points)))
315 data.shape_values.reinit(this->n_dofs_per_cell(), n_q_points);
316
317 if (update_flags & update_gradients)
318 data.shape_gradients.reinit(this->n_dofs_per_cell(), n_q_points);
319
320 if (update_flags & update_hessians)
321 data.shape_hessians.reinit(this->n_dofs_per_cell(), n_q_points);
322
323 if (update_flags & update_3rd_derivatives)
324 data.shape_3rd_derivatives.reinit(this->n_dofs_per_cell(), n_q_points);
325
326 // next already fill those fields of which we have information by
327 // now. note that the shape gradients are only those on the unit
328 // cell, and need to be transformed when visiting an actual cell
329 if (update_flags & (update_values | update_gradients | update_hessians |
331 for (unsigned int i = 0; i < n_q_points; ++i)
332 {
333 poly_space->evaluate(quadrature.point(i),
334 values,
335 grads,
336 grad_grads,
337 third_derivatives,
338 fourth_derivatives);
339
340 // the values of shape functions at quadrature points don't change.
341 // consequently, write these values right into the output array if
342 // we can, i.e., if the output array has the correct size. this is
343 // the case on cells. on faces, we already precompute data on *all*
344 // faces and subfaces, but we later on copy only a portion of it
345 // into the output object; in that case, copy the data from all
346 // faces into the scratch object
347 if (update_flags & update_values)
348 if (output_data.shape_values.n_rows() > 0)
349 {
350 if (output_data.shape_values.n_cols() == n_q_points)
351 for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
352 output_data.shape_values[k][i] = values[k];
353 else
354 for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
355 data.shape_values[k][i] = values[k];
356 }
357
358 // for everything else, derivatives need to be transformed,
359 // so we write them into our scratch space and only later
360 // copy stuff into where FEValues wants it
361 if (update_flags & update_gradients)
362 for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
363 data.shape_gradients[k][i] = grads[k];
364
365 if (update_flags & update_hessians)
366 for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
367 data.shape_hessians[k][i] = grad_grads[k];
368
369 if (update_flags & update_3rd_derivatives)
370 for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
371 data.shape_3rd_derivatives[k][i] = third_derivatives[k];
372 }
373 return data_ptr;
374 }
375
376 virtual void
379 const CellSimilarity::Similarity cell_similarity,
380 const Quadrature<dim> & quadrature,
381 const Mapping<dim, spacedim> & mapping,
382 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
384 & mapping_data,
385 const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
387 spacedim>
388 &output_data) const override;
389
390 using FiniteElement<dim, spacedim>::fill_fe_face_values;
391
392 virtual void
395 const unsigned int face_no,
396 const hp::QCollection<dim - 1> & quadrature,
397 const Mapping<dim, spacedim> & mapping,
398 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
400 & mapping_data,
401 const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
403 spacedim>
404 &output_data) const override;
405
406 virtual void
409 const unsigned int face_no,
410 const unsigned int sub_no,
411 const Quadrature<dim - 1> & quadrature,
412 const Mapping<dim, spacedim> & mapping,
413 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
415 & mapping_data,
416 const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
418 spacedim>
419 &output_data) const override;
420
427 class InternalData : public FiniteElement<dim, spacedim>::InternalDataBase
428 {
429 public:
440
451
462
473 };
474
491 void
494 &output_data,
496 & mapping_data,
497 const unsigned int n_q_points) const;
498
521 void
524 &output_data,
526 & mapping_data,
527 const unsigned int n_q_points) const;
528
529
533 const std::unique_ptr<ScalarPolynomialsBase<dim>> poly_space;
534};
535
539
540#endif
Table< 2, double > shape_values
Definition fe_poly.h:439
Table< 2, Tensor< 3, dim > > shape_3rd_derivatives
Definition fe_poly.h:472
Table< 2, Tensor< 2, dim > > shape_hessians
Definition fe_poly.h:461
Table< 2, Tensor< 1, dim > > shape_gradients
Definition fe_poly.h:450
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const override
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
FE_Poly(const FE_Poly &fe)
virtual Tensor< 3, dim > shape_3rd_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
std::vector< unsigned int > get_poly_space_numbering_inverse() const
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) 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
Definition fe_poly.h:263
FE_Poly(const ScalarPolynomialsBase< dim > &poly_space, const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags, const std::vector< ComponentMask > &nonzero_components)
void correct_hessians(internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const unsigned int n_q_points) const
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 Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const override
unsigned int get_degree() const
const ScalarPolynomialsBase< dim > & get_poly_space() const
virtual Tensor< 4, dim > shape_4th_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
void correct_third_derivatives(internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const unsigned int n_q_points) const
virtual double shape_value(const unsigned int i, const Point< dim > &p) const override
const std::unique_ptr< ScalarPolynomialsBase< dim > > poly_space
Definition fe_poly.h:533
virtual Tensor< 3, dim > shape_3rd_derivative(const unsigned int i, const Point< dim > &p) 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
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
virtual Tensor< 4, dim > shape_4th_derivative(const unsigned int i, const Point< dim > &p) const override
std::vector< unsigned int > get_poly_space_numbering() const
virtual std::size_t memory_consumption() 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 Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
unsigned int n_dofs_per_cell() const
const std::vector< bool > restriction_is_additive_flags
Definition fe.h:2570
const std::vector< ComponentMask > nonzero_components
Definition fe.h:2579
Abstract base class for mapping classes.
Definition mapping.h:317
Definition point.h:112
const Point< dim > & point(const unsigned int i) const
unsigned int size() const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
UpdateFlags
@ update_hessians
Second derivatives of shape functions.
@ update_values
Shape function values.
@ update_3rd_derivatives
Third derivatives of shape functions.
@ update_gradients
Shape function gradients.