Reference documentation for deal.II version Git d3aed38b93 2021-10-28 13:33:27 +0200
\(\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_poly.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2004 - 2021 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 
33 
73 template <int dim, int spacedim = dim>
74 class FE_Poly : public FiniteElement<dim, spacedim>
75 {
76 public:
81  const FiniteElementData<dim> & fe_data,
82  const std::vector<bool> & restriction_is_additive_flags,
83  const std::vector<ComponentMask> &nonzero_components);
84 
88  FE_Poly(const FE_Poly &fe);
89 
94  unsigned int
95  get_degree() const;
96 
97  // for documentation, see the FiniteElement base class
98  virtual UpdateFlags
99  requires_update_flags(const UpdateFlags update_flags) const override;
100 
105  get_poly_space() const;
106 
116  std::vector<unsigned int>
117  get_poly_space_numbering() const;
118 
125  std::vector<unsigned int>
127 
133  virtual double
134  shape_value(const unsigned int i, const Point<dim> &p) const override;
135 
146  virtual double
147  shape_value_component(const unsigned int i,
148  const Point<dim> & p,
149  const unsigned int component) const override;
150 
156  virtual Tensor<1, dim>
157  shape_grad(const unsigned int i, const Point<dim> &p) const override;
158 
169  virtual Tensor<1, dim>
170  shape_grad_component(const unsigned int i,
171  const Point<dim> & p,
172  const unsigned int component) const override;
173 
179  virtual Tensor<2, dim>
180  shape_grad_grad(const unsigned int i, const Point<dim> &p) const override;
181 
192  virtual Tensor<2, dim>
193  shape_grad_grad_component(const unsigned int i,
194  const Point<dim> & p,
195  const unsigned int component) const override;
196 
202  virtual Tensor<3, dim>
203  shape_3rd_derivative(const unsigned int i,
204  const Point<dim> & p) const override;
205 
216  virtual Tensor<3, dim>
217  shape_3rd_derivative_component(const unsigned int i,
218  const Point<dim> & p,
219  const unsigned int component) const override;
220 
226  virtual Tensor<4, dim>
227  shape_4th_derivative(const unsigned int i,
228  const Point<dim> & p) const override;
229 
240  virtual Tensor<4, dim>
241  shape_4th_derivative_component(const unsigned int i,
242  const Point<dim> & p,
243  const unsigned int component) const override;
244 
248  virtual std::size_t
249  memory_consumption() const override;
250 
251 protected:
252  /*
253  * NOTE: The following function has its definition inlined into the class
254  * declaration because we otherwise run into a compiler error with MS Visual
255  * Studio.
256  */
257 
258 
259  virtual std::unique_ptr<
262  const UpdateFlags update_flags,
263  const Mapping<dim, spacedim> &mapping,
264  const Quadrature<dim> & quadrature,
266  spacedim>
267  &output_data) const override
268  {
269  (void)mapping;
270 
271  // generate a new data object and
272  // initialize some fields
273  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
274  data_ptr = std::make_unique<InternalData>();
275  auto &data = dynamic_cast<InternalData &>(*data_ptr);
276  data.update_each = requires_update_flags(update_flags);
277 
278  const unsigned int n_q_points = quadrature.size();
279 
280  // initialize some scratch arrays. we need them for the underlying
281  // polynomial to put the values and derivatives of shape functions
282  // to put there, depending on what the user requested
283  std::vector<double> values(
284  update_flags & update_values ? this->n_dofs_per_cell() : 0);
285  std::vector<Tensor<1, dim>> grads(
286  update_flags & update_gradients ? this->n_dofs_per_cell() : 0);
287  std::vector<Tensor<2, dim>> grad_grads(
288  update_flags & update_hessians ? this->n_dofs_per_cell() : 0);
289  std::vector<Tensor<3, dim>> third_derivatives(
290  update_flags & update_3rd_derivatives ? this->n_dofs_per_cell() : 0);
291  std::vector<Tensor<4, dim>>
292  fourth_derivatives; // won't be needed, so leave empty
293 
294  // now also initialize fields the fields of this class's own
295  // temporary storage, depending on what we need for the given
296  // update flags.
297  //
298  // there is one exception from the rule: if we are dealing with
299  // cells (i.e., if this function is not called via
300  // get_(sub)face_data()), then we can already store things in the
301  // final location where FEValues::reinit() later wants to see
302  // things. we then don't need the intermediate space. we determine
303  // whether we are on a cell by asking whether the number of
304  // elements in the output array equals the number of quadrature
305  // points (yes, it's a cell) or not (because in that case the
306  // number of quadrature points we use here equals the number of
307  // quadrature points summed over *all* faces or subfaces, whereas
308  // the number of output slots equals the number of quadrature
309  // points on only *one* face)
310  if ((update_flags & update_values) &&
311  !((output_data.shape_values.n_rows() > 0) &&
312  (output_data.shape_values.n_cols() == n_q_points)))
313  data.shape_values.reinit(this->n_dofs_per_cell(), n_q_points);
314 
315  if (update_flags & update_gradients)
316  data.shape_gradients.reinit(this->n_dofs_per_cell(), n_q_points);
317 
318  if (update_flags & update_hessians)
319  data.shape_hessians.reinit(this->n_dofs_per_cell(), n_q_points);
320 
321  if (update_flags & update_3rd_derivatives)
322  data.shape_3rd_derivatives.reinit(this->n_dofs_per_cell(), n_q_points);
323 
324  // next already fill those fields of which we have information by
325  // now. note that the shape gradients are only those on the unit
326  // cell, and need to be transformed when visiting an actual cell
327  if (update_flags & (update_values | update_gradients | update_hessians |
328  update_3rd_derivatives))
329  for (unsigned int i = 0; i < n_q_points; ++i)
330  {
331  poly_space->evaluate(quadrature.point(i),
332  values,
333  grads,
334  grad_grads,
335  third_derivatives,
336  fourth_derivatives);
337 
338  // the values of shape functions at quadrature points don't change.
339  // consequently, write these values right into the output array if
340  // we can, i.e., if the output array has the correct size. this is
341  // the case on cells. on faces, we already precompute data on *all*
342  // faces and subfaces, but we later on copy only a portion of it
343  // into the output object; in that case, copy the data from all
344  // faces into the scratch object
345  if (update_flags & update_values)
346  if (output_data.shape_values.n_rows() > 0)
347  {
348  if (output_data.shape_values.n_cols() == n_q_points)
349  for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
350  output_data.shape_values[k][i] = values[k];
351  else
352  for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
353  data.shape_values[k][i] = values[k];
354  }
355 
356  // for everything else, derivatives need to be transformed,
357  // so we write them into our scratch space and only later
358  // copy stuff into where FEValues wants it
359  if (update_flags & update_gradients)
360  for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
361  data.shape_gradients[k][i] = grads[k];
362 
363  if (update_flags & update_hessians)
364  for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
365  data.shape_hessians[k][i] = grad_grads[k];
366 
367  if (update_flags & update_3rd_derivatives)
368  for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
369  data.shape_3rd_derivatives[k][i] = third_derivatives[k];
370  }
371  return data_ptr;
372  }
373 
374  virtual void
376  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
377  const CellSimilarity::Similarity cell_similarity,
378  const Quadrature<dim> & quadrature,
379  const Mapping<dim, spacedim> & mapping,
380  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
381  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
382  spacedim>
383  & mapping_data,
384  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
386  spacedim>
387  &output_data) const override;
388 
390 
391  virtual void
393  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
394  const unsigned int face_no,
395  const hp::QCollection<dim - 1> & quadrature,
396  const Mapping<dim, spacedim> & mapping,
397  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
398  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
399  spacedim>
400  & mapping_data,
401  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
403  spacedim>
404  &output_data) const override;
405 
406  virtual void
408  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
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,
414  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
415  spacedim>
416  & mapping_data,
417  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
419  spacedim>
420  &output_data) const override;
421 
428  class InternalData : public FiniteElement<dim, spacedim>::InternalDataBase
429  {
430  public:
441 
452 
463 
474  };
475 
492  void
495  &output_data,
497  & mapping_data,
498  const unsigned int n_q_points) const;
499 
522  void
525  &output_data,
527  & mapping_data,
528  const unsigned int n_q_points) const;
529 
530 
534  const std::unique_ptr<ScalarPolynomialsBase<dim>> poly_space;
535 };
536 
540 
541 #endif
Shape function values.
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)
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const override
virtual std::size_t memory_consumption() const override
const Point< dim > & point(const unsigned int i) const
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const override
const ScalarPolynomialsBase< dim > & get_poly_space() const
Table< 2, double > shape_values
Definition: fe_poly.h:440
virtual double shape_value(const unsigned int i, const Point< dim > &p) 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
Third derivatives of shape functions.
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) 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
UpdateFlags
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
Abstract base class for mapping classes.
Definition: mapping.h:303
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
virtual Tensor< 3, dim > shape_3rd_derivative(const unsigned int i, const Point< dim > &p) const override
Second derivatives of shape functions.
unsigned int size() const
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
virtual Tensor< 4, dim > shape_4th_derivative(const unsigned int i, const Point< dim > &p) const override
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_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
unsigned int get_degree() const
unsigned int n_dofs_per_cell() const
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:261
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
Shape function gradients.
virtual Tensor< 4, dim > shape_4th_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
const std::vector< ComponentMask > nonzero_components
Definition: fe.h:2581
Table< 2, Tensor< 2, dim > > shape_hessians
Definition: fe_poly.h:462
Table< 2, Tensor< 1, dim > > shape_gradients
Definition: fe_poly.h:451
Table< 2, Tensor< 3, dim > > shape_3rd_derivatives
Definition: fe_poly.h:473
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
std::vector< unsigned int > get_poly_space_numbering() const
const std::unique_ptr< ScalarPolynomialsBase< dim > > poly_space
Definition: fe_poly.h:534
const std::vector< bool > restriction_is_additive_flags
Definition: fe.h:2572