Reference documentation for deal.II version Git 840703ecf3 2019-10-14 17:15:20 -0400
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
fe_poly.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2004 - 2019 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/quadrature.h>
21 #include <deal.II/base/std_cxx14/memory.h>
22 
23 #include <deal.II/fe/fe.h>
24 
25 DEAL_II_NAMESPACE_OPEN
26 
29 
71 template <class PolynomialType,
72  int dim = PolynomialType::dimension,
73  int spacedim = dim>
74 class FE_Poly : public FiniteElement<dim, spacedim>
75 {
76 public:
80  FE_Poly(const PolynomialType & poly_space,
81  const FiniteElementData<dim> & fe_data,
82  const std::vector<bool> & restriction_is_additive_flags,
83  const std::vector<ComponentMask> &nonzero_components);
84 
89  unsigned int
90  get_degree() const;
91 
92  // for documentation, see the FiniteElement base class
93  virtual UpdateFlags
94  requires_update_flags(const UpdateFlags update_flags) const override;
95 
101  std::vector<unsigned int>
102  get_poly_space_numbering() const;
103 
108  std::vector<unsigned int>
110 
116  virtual double
117  shape_value(const unsigned int i, const Point<dim> &p) const override;
118 
129  virtual double
130  shape_value_component(const unsigned int i,
131  const Point<dim> & p,
132  const unsigned int component) const override;
133 
139  virtual Tensor<1, dim>
140  shape_grad(const unsigned int i, const Point<dim> &p) const override;
141 
152  virtual Tensor<1, dim>
153  shape_grad_component(const unsigned int i,
154  const Point<dim> & p,
155  const unsigned int component) const override;
156 
162  virtual Tensor<2, dim>
163  shape_grad_grad(const unsigned int i, const Point<dim> &p) const override;
164 
175  virtual Tensor<2, dim>
176  shape_grad_grad_component(const unsigned int i,
177  const Point<dim> & p,
178  const unsigned int component) const override;
179 
185  virtual Tensor<3, dim>
186  shape_3rd_derivative(const unsigned int i,
187  const Point<dim> & p) const override;
188 
199  virtual Tensor<3, dim>
200  shape_3rd_derivative_component(const unsigned int i,
201  const Point<dim> & p,
202  const unsigned int component) const override;
203 
209  virtual Tensor<4, dim>
210  shape_4th_derivative(const unsigned int i,
211  const Point<dim> & p) const override;
212 
223  virtual Tensor<4, dim>
224  shape_4th_derivative_component(const unsigned int i,
225  const Point<dim> & p,
226  const unsigned int component) const override;
227 
228 protected:
229  /*
230  * NOTE: The following function has its definition inlined into the class
231  * declaration because we otherwise run into a compiler error with MS Visual
232  * Studio.
233  */
234 
235 
236  virtual std::unique_ptr<
239  const UpdateFlags update_flags,
240  const Mapping<dim, spacedim> & /*mapping*/,
241  const Quadrature<dim> &quadrature,
243  spacedim>
244  &output_data) const override
245  {
246  // generate a new data object and
247  // initialize some fields
248  std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
249  data_ptr = std_cxx14::make_unique<InternalData>();
250  auto &data = dynamic_cast<InternalData &>(*data_ptr);
251  data.update_each = requires_update_flags(update_flags);
252 
253  const unsigned int n_q_points = quadrature.size();
254 
255  // initialize some scratch arrays. we need them for the underlying
256  // polynomial to put the values and derivatives of shape functions
257  // to put there, depending on what the user requested
258  std::vector<double> values(
259  update_flags & update_values ? this->dofs_per_cell : 0);
260  std::vector<Tensor<1, dim>> grads(
261  update_flags & update_gradients ? this->dofs_per_cell : 0);
262  std::vector<Tensor<2, dim>> grad_grads(
263  update_flags & update_hessians ? this->dofs_per_cell : 0);
264  std::vector<Tensor<3, dim>> third_derivatives(
265  update_flags & update_3rd_derivatives ? this->dofs_per_cell : 0);
266  std::vector<Tensor<4, dim>>
267  fourth_derivatives; // won't be needed, so leave empty
268 
269  // now also initialize fields the fields of this class's own
270  // temporary storage, depending on what we need for the given
271  // update flags.
272  //
273  // there is one exception from the rule: if we are dealing with
274  // cells (i.e., if this function is not called via
275  // get_(sub)face_data()), then we can already store things in the
276  // final location where FEValues::reinit() later wants to see
277  // things. we then don't need the intermediate space. we determine
278  // whether we are on a cell by asking whether the number of
279  // elements in the output array equals the number of quadrature
280  // points (yes, it's a cell) or not (because in that case the
281  // number of quadrature points we use here equals the number of
282  // quadrature points summed over *all* faces or subfaces, whereas
283  // the number of output slots equals the number of quadrature
284  // points on only *one* face)
285  if ((update_flags & update_values) &&
286  !((output_data.shape_values.n_rows() > 0) &&
287  (output_data.shape_values.n_cols() == n_q_points)))
288  data.shape_values.reinit(this->dofs_per_cell, n_q_points);
289 
290  if (update_flags & update_gradients)
291  data.shape_gradients.reinit(this->dofs_per_cell, n_q_points);
292 
293  if (update_flags & update_hessians)
294  data.shape_hessians.reinit(this->dofs_per_cell, n_q_points);
295 
296  if (update_flags & update_3rd_derivatives)
297  data.shape_3rd_derivatives.reinit(this->dofs_per_cell, n_q_points);
298 
299  // next already fill those fields of which we have information by
300  // now. note that the shape gradients are only those on the unit
301  // cell, and need to be transformed when visiting an actual cell
302  if (update_flags & (update_values | update_gradients | update_hessians |
303  update_3rd_derivatives))
304  for (unsigned int i = 0; i < n_q_points; ++i)
305  {
306  poly_space.evaluate(quadrature.point(i),
307  values,
308  grads,
309  grad_grads,
310  third_derivatives,
311  fourth_derivatives);
312 
313  // the values of shape functions at quadrature points don't change.
314  // consequently, write these values right into the output array if
315  // we can, i.e., if the output array has the correct size. this is
316  // the case on cells. on faces, we already precompute data on *all*
317  // faces and subfaces, but we later on copy only a portion of it
318  // into the output object; in that case, copy the data from all
319  // faces into the scratch object
320  if (update_flags & update_values)
321  if (output_data.shape_values.n_rows() > 0)
322  {
323  if (output_data.shape_values.n_cols() == n_q_points)
324  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
325  output_data.shape_values[k][i] = values[k];
326  else
327  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
328  data.shape_values[k][i] = values[k];
329  }
330 
331  // for everything else, derivatives need to be transformed,
332  // so we write them into our scratch space and only later
333  // copy stuff into where FEValues wants it
334  if (update_flags & update_gradients)
335  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
336  data.shape_gradients[k][i] = grads[k];
337 
338  if (update_flags & update_hessians)
339  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
340  data.shape_hessians[k][i] = grad_grads[k];
341 
342  if (update_flags & update_3rd_derivatives)
343  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
344  data.shape_3rd_derivatives[k][i] = third_derivatives[k];
345  }
346  return data_ptr;
347  }
348 
349  virtual void
350  fill_fe_values(
351  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
352  const CellSimilarity::Similarity cell_similarity,
353  const Quadrature<dim> & quadrature,
354  const Mapping<dim, spacedim> & mapping,
355  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
356  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
357  spacedim>
358  & mapping_data,
359  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
361  spacedim>
362  &output_data) const override;
363 
364  virtual void
365  fill_fe_face_values(
366  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
367  const unsigned int face_no,
368  const Quadrature<dim - 1> & quadrature,
369  const Mapping<dim, spacedim> & mapping,
370  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
371  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
372  spacedim>
373  & mapping_data,
374  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
376  spacedim>
377  &output_data) const override;
378 
379  virtual void
380  fill_fe_subface_values(
381  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
382  const unsigned int face_no,
383  const unsigned int sub_no,
384  const Quadrature<dim - 1> & quadrature,
385  const Mapping<dim, spacedim> & mapping,
386  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
387  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
388  spacedim>
389  & mapping_data,
390  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
392  spacedim>
393  &output_data) const override;
394 
401  class InternalData : public FiniteElement<dim, spacedim>::InternalDataBase
402  {
403  public:
414 
425 
436 
447  };
448 
471  void
474  &output_data,
476  & mapping_data,
477  const unsigned int n_q_points,
478  const unsigned int dof) const;
479 
484  PolynomialType poly_space;
485 };
486 
489 DEAL_II_NAMESPACE_CLOSE
490 
491 #endif
Shape function values.
virtual double shape_value(const unsigned int i, const Point< dim > &p) const override
Table< 2, Tensor< 2, dim > > shape_hessians
Definition: fe_poly.h:435
virtual Tensor< 4, dim > shape_4th_derivative(const unsigned int i, const Point< dim > &p) const override
const Point< dim > & point(const unsigned int i) const
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const override
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
Table< 2, Tensor< 1, dim > > shape_gradients
Definition: fe_poly.h:424
FE_Poly(const PolynomialType &poly_space, const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags, const std::vector< ComponentMask > &nonzero_components)
virtual std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &, const Quadrature< dim > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
Definition: fe_poly.h:238
Third derivatives of shape functions.
UpdateFlags
Abstract base class for mapping classes.
Definition: mapping.h:302
Table< 2, Tensor< 3, dim > > shape_3rd_derivatives
Definition: fe_poly.h:446
Second derivatives of shape functions.
unsigned int size() const
const unsigned int dofs_per_cell
Definition: fe_base.h:282
virtual Tensor< 3, dim > shape_3rd_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
Table< 2, double > shape_values
Definition: fe_poly.h:413
unsigned int get_degree() const
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 unsigned int dof) const
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
virtual Tensor< 4, dim > shape_4th_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
Shape function gradients.
virtual Tensor< 3, dim > shape_3rd_derivative(const unsigned int i, const Point< dim > &p) const override
const std::vector< ComponentMask > nonzero_components
Definition: fe.h:2618
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
PolynomialType poly_space
Definition: fe_poly.h:484
std::vector< unsigned int > get_poly_space_numbering_inverse() const
std::vector< unsigned int > get_poly_space_numbering() const
const std::vector< bool > restriction_is_additive_flags
Definition: fe.h:2609
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override