Reference documentation for deal.II version Git 6cf3cbfd16 2020-07-12 08:18:15 -0400
\(\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_tensor.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2005 - 2020 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_tensor_h
17 #define dealii_fe_poly_tensor_h
18 
19 
20 #include <deal.II/base/config.h>
21 
26 
27 #include <deal.II/fe/fe.h>
28 
30 
31 #include <memory>
32 
34 
140 template <int dim, int spacedim = dim>
141 class FE_PolyTensor : public FiniteElement<dim, spacedim>
142 {
143 public:
150  FE_PolyTensor(const TensorPolynomialsBase<dim> &polynomials,
151  const FiniteElementData<dim> & fe_data,
152  const std::vector<bool> & restriction_is_additive_flags,
153  const std::vector<ComponentMask> &nonzero_components);
154 
155 
159  FE_PolyTensor(const FE_PolyTensor &fe);
160 
161  // for documentation, see the FiniteElement base class
162  virtual UpdateFlags
163  requires_update_flags(const UpdateFlags update_flags) const override;
164 
171  virtual double
172  shape_value(const unsigned int i, const Point<dim> &p) const override;
173 
174  // documentation inherited from the base class
175  virtual double
176  shape_value_component(const unsigned int i,
177  const Point<dim> & p,
178  const unsigned int component) const override;
179 
186  virtual Tensor<1, dim>
187  shape_grad(const unsigned int i, const Point<dim> &p) const override;
188 
189  // documentation inherited from the base class
190  virtual Tensor<1, dim>
191  shape_grad_component(const unsigned int i,
192  const Point<dim> & p,
193  const unsigned int component) const override;
194 
201  virtual Tensor<2, dim>
202  shape_grad_grad(const unsigned int i, const Point<dim> &p) const override;
203 
204  // documentation inherited from the base class
205  virtual Tensor<2, dim>
206  shape_grad_grad_component(const unsigned int i,
207  const Point<dim> & p,
208  const unsigned int component) const override;
209 
210 protected:
218  std::vector<MappingKind> mapping_kind;
219 
224  bool
225  single_mapping_kind() const;
226 
231  get_mapping_kind(const unsigned int i) const;
232 
233  /* NOTE: The following function has its definition inlined into the class
234  declaration because we otherwise run into a compiler error with MS Visual
235  Studio. */
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::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  // some scratch arrays
256  std::vector<Tensor<1, dim>> values(0);
257  std::vector<Tensor<2, dim>> grads(0);
258  std::vector<Tensor<3, dim>> grad_grads(0);
259  std::vector<Tensor<4, dim>> third_derivatives(0);
260  std::vector<Tensor<5, dim>> fourth_derivatives(0);
261 
262  if (update_flags & (update_values | update_gradients | update_hessians))
263  data.sign_change.resize(this->dofs_per_cell);
264 
265  // initialize fields only if really
266  // necessary. otherwise, don't
267  // allocate memory
268 
269  const bool update_transformed_shape_values =
270  std::any_of(this->mapping_kind.begin(),
271  this->mapping_kind.end(),
272  [](const MappingKind t) { return t != mapping_none; });
273 
274  const bool update_transformed_shape_grads =
275  std::any_of(this->mapping_kind.begin(),
276  this->mapping_kind.end(),
277  [](const MappingKind t) {
278  return (t == mapping_raviart_thomas || t == mapping_piola ||
280  });
281 
282  const bool update_transformed_shape_hessian_tensors =
283  update_transformed_shape_values;
284 
285  if (update_flags & update_values)
286  {
287  values.resize(this->dofs_per_cell);
288  data.shape_values.reinit(this->dofs_per_cell, n_q_points);
289  if (update_transformed_shape_values)
290  data.transformed_shape_values.resize(n_q_points);
291  }
292 
293  if (update_flags & update_gradients)
294  {
295  grads.resize(this->dofs_per_cell);
296  data.shape_grads.reinit(this->dofs_per_cell, n_q_points);
297  data.transformed_shape_grads.resize(n_q_points);
298 
299  if (update_transformed_shape_grads)
300  data.untransformed_shape_grads.resize(n_q_points);
301  }
302 
303  if (update_flags & update_hessians)
304  {
305  grad_grads.resize(this->dofs_per_cell);
306  data.shape_grad_grads.reinit(this->dofs_per_cell, n_q_points);
307  data.transformed_shape_hessians.resize(n_q_points);
308  if (update_transformed_shape_hessian_tensors)
309  data.untransformed_shape_hessian_tensors.resize(n_q_points);
310  }
311 
312  // Compute shape function values
313  // and derivatives and hessians on
314  // the reference cell.
315  // Make sure, that for the
316  // node values N_i holds
317  // N_i(v_j)=\delta_ij for all basis
318  // functions v_j
319  if (update_flags & (update_values | update_gradients))
320  for (unsigned int k = 0; k < n_q_points; ++k)
321  {
322  poly_space->evaluate(quadrature.point(k),
323  values,
324  grads,
325  grad_grads,
326  third_derivatives,
327  fourth_derivatives);
328 
329  if (update_flags & update_values)
330  {
331  if (inverse_node_matrix.n_cols() == 0)
332  for (unsigned int i = 0; i < this->dofs_per_cell; ++i)
333  data.shape_values[i][k] = values[i];
334  else
335  for (unsigned int i = 0; i < this->dofs_per_cell; ++i)
336  {
337  Tensor<1, dim> add_values;
338  for (unsigned int j = 0; j < this->dofs_per_cell; ++j)
339  add_values += inverse_node_matrix(j, i) * values[j];
340  data.shape_values[i][k] = add_values;
341  }
342  }
343 
344  if (update_flags & update_gradients)
345  {
346  if (inverse_node_matrix.n_cols() == 0)
347  for (unsigned int i = 0; i < this->dofs_per_cell; ++i)
348  data.shape_grads[i][k] = grads[i];
349  else
350  for (unsigned int i = 0; i < this->dofs_per_cell; ++i)
351  {
352  Tensor<2, dim> add_grads;
353  for (unsigned int j = 0; j < this->dofs_per_cell; ++j)
354  add_grads += inverse_node_matrix(j, i) * grads[j];
355  data.shape_grads[i][k] = add_grads;
356  }
357  }
358 
359  if (update_flags & update_hessians)
360  {
361  if (inverse_node_matrix.n_cols() == 0)
362  for (unsigned int i = 0; i < this->dofs_per_cell; ++i)
363  data.shape_grad_grads[i][k] = grad_grads[i];
364  else
365  for (unsigned int i = 0; i < this->dofs_per_cell; ++i)
366  {
367  Tensor<3, dim> add_grad_grads;
368  for (unsigned int j = 0; j < this->dofs_per_cell; ++j)
369  add_grad_grads +=
370  inverse_node_matrix(j, i) * grad_grads[j];
371  data.shape_grad_grads[i][k] = add_grad_grads;
372  }
373  }
374  }
375  return data_ptr;
376  }
377 
378  virtual void
380  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
381  const CellSimilarity::Similarity cell_similarity,
382  const Quadrature<dim> & quadrature,
383  const Mapping<dim, spacedim> & mapping,
384  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
385  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
386  spacedim>
387  & mapping_data,
388  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
390  spacedim>
391  &output_data) const override;
392 
393  virtual void
395  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
396  const unsigned int face_no,
397  const Quadrature<dim - 1> & quadrature,
398  const Mapping<dim, spacedim> & mapping,
399  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
400  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
401  spacedim>
402  & mapping_data,
403  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
405  spacedim>
406  &output_data) const override;
407 
408  virtual void
410  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
411  const unsigned int face_no,
412  const unsigned int sub_no,
413  const Quadrature<dim - 1> & quadrature,
414  const Mapping<dim, spacedim> & mapping,
415  const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
416  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
417  spacedim>
418  & mapping_data,
419  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
421  spacedim>
422  &output_data) const override;
423 
433  class InternalData : public FiniteElement<dim, spacedim>::InternalDataBase
434  {
435  public:
441 
448 
455 
459  mutable std::vector<double> sign_change;
460  mutable std::vector<Tensor<1, spacedim>> transformed_shape_values;
461  // for shape_gradient computations
462  mutable std::vector<Tensor<2, spacedim>> transformed_shape_grads;
463  mutable std::vector<Tensor<2, dim>> untransformed_shape_grads;
464  // for shape_hessian computations
465  mutable std::vector<Tensor<3, spacedim>> transformed_shape_hessians;
466  mutable std::vector<Tensor<3, dim>> untransformed_shape_hessian_tensors;
467  };
468 
469 
470 
475  const std::unique_ptr<const TensorPolynomialsBase<dim>> poly_space;
476 
489 
493  mutable std::mutex cache_mutex;
494 
501 
505  mutable std::vector<Tensor<1, dim>> cached_values;
506 
510  mutable std::vector<Tensor<2, dim>> cached_grads;
511 
516  mutable std::vector<Tensor<3, dim>> cached_grad_grads;
517 };
518 
520 
521 #endif
Shape function values.
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
std::vector< Tensor< 2, dim > > cached_grads
MappingKind get_mapping_kind(const unsigned int i) const
bool single_mapping_kind() const
std::vector< Tensor< 1, dim > > cached_values
std::vector< double > sign_change
const Point< dim > & point(const unsigned int i) const
FullMatrix< double > inverse_node_matrix
Table< 2, DerivativeForm< 1, dim, spacedim > > shape_grads
MappingKind
Definition: mapping.h:62
std::vector< Tensor< 2, dim > > untransformed_shape_grads
std::vector< Tensor< 3, spacedim > > transformed_shape_hessians
Table< 2, DerivativeForm< 2, dim, spacedim > > shape_grad_grads
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< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
std::vector< Tensor< 3, dim > > untransformed_shape_hessian_tensors
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 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 > &) const override
Abstract base class for mapping classes.
Definition: mapping.h:301
std::mutex cache_mutex
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
Point< dim > cached_point
std::vector< Tensor< 1, spacedim > > transformed_shape_values
std::vector< Tensor< 3, dim > > cached_grad_grads
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const override
virtual double shape_value(const unsigned int i, const Point< dim > &p) const override
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const override
Second derivatives of shape functions.
unsigned int size() const
const unsigned int dofs_per_cell
Definition: fe_base.h:280
const std::unique_ptr< const TensorPolynomialsBase< dim > > poly_space
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:362
Shape function gradients.
Table< 2, Tensor< 1, dim > > shape_values
const std::vector< ComponentMask > nonzero_components
Definition: fe.h:2570
Definition: table.h:687
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_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
FE_PolyTensor(const TensorPolynomialsBase< dim > &polynomials, const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags, const std::vector< ComponentMask > &nonzero_components)
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
std::vector< Tensor< 2, spacedim > > transformed_shape_grads
const std::vector< bool > restriction_is_additive_flags
Definition: fe.h:2561
std::vector< MappingKind > mapping_kind