Reference documentation for deal.II version 9.2.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\}}\)
fe_poly.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2004 - 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 
24 
25 #include <deal.II/fe/fe_poly.h>
26 #include <deal.II/fe/fe_poly.templates.h>
27 #include <deal.II/fe/fe_values.h>
28 
30 
31 
32 template <>
33 void
34 FE_Poly<TensorProductPolynomials<1>, 1, 2>::fill_fe_values(
36  const CellSimilarity::Similarity cell_similarity,
37  const Quadrature<1> & quadrature,
38  const Mapping<1, 2> & mapping,
39  const Mapping<1, 2>::InternalDataBase &mapping_internal,
40  const ::internal::FEValuesImplementation::MappingRelatedData<1, 2>
41  & mapping_data,
42  const FiniteElement<1, 2>::InternalDataBase &fe_internal,
44  &output_data) const
45 {
46  // convert data object to internal data for this class. fails with an
47  // exception if that is not possible
48  Assert(dynamic_cast<const InternalData *>(&fe_internal) != nullptr,
50  const InternalData &fe_data =
51  static_cast<const InternalData &>(fe_internal); // NOLINT
52 
53  const bool need_to_correct_higher_derivatives =
54  higher_derivatives_need_correcting(mapping,
55  mapping_data,
56  quadrature.size(),
57  fe_data.update_each);
58 
59  // transform gradients and higher derivatives. there is nothing to do
60  // for values since we already emplaced them into output_data when
61  // we were in get_data()
62  if (fe_data.update_each & update_gradients &&
63  cell_similarity != CellSimilarity::translation)
64  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
65  mapping.transform(make_array_view(fe_data.shape_gradients, k),
67  mapping_internal,
68  make_array_view(output_data.shape_gradients, k));
69 
70  if (fe_data.update_each & update_hessians &&
71  cell_similarity != CellSimilarity::translation)
72  {
73  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
74  mapping.transform(make_array_view(fe_data.shape_hessians, k),
76  mapping_internal,
77  make_array_view(output_data.shape_hessians, k));
78 
79  if (need_to_correct_higher_derivatives)
80  correct_hessians(output_data, mapping_data, quadrature.size());
81  }
82 
83  if (fe_data.update_each & update_3rd_derivatives &&
84  cell_similarity != CellSimilarity::translation)
85  {
86  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
87  mapping.transform(make_array_view(fe_data.shape_3rd_derivatives, k),
89  mapping_internal,
91  k));
92 
93  if (need_to_correct_higher_derivatives)
94  correct_third_derivatives(output_data, mapping_data, quadrature.size());
95  }
96 }
97 
98 
99 
100 template <>
101 void
102 FE_Poly<TensorProductPolynomials<2>, 2, 3>::fill_fe_values(
104  const CellSimilarity::Similarity cell_similarity,
105  const Quadrature<2> & quadrature,
106  const Mapping<2, 3> & mapping,
107  const Mapping<2, 3>::InternalDataBase &mapping_internal,
108  const ::internal::FEValuesImplementation::MappingRelatedData<2, 3>
109  & mapping_data,
110  const FiniteElement<2, 3>::InternalDataBase &fe_internal,
112  &output_data) const
113 {
114  // convert data object to internal data for this class. fails with an
115  // exception if that is not possible
116  Assert(dynamic_cast<const InternalData *>(&fe_internal) != nullptr,
117  ExcInternalError());
118  const InternalData &fe_data =
119  static_cast<const InternalData &>(fe_internal); // NOLINT
120 
121  const bool need_to_correct_higher_derivatives =
122  higher_derivatives_need_correcting(mapping,
123  mapping_data,
124  quadrature.size(),
125  fe_data.update_each);
126 
127  // transform gradients and higher derivatives. there is nothing to do
128  // for values since we already emplaced them into output_data when
129  // we were in get_data()
130  if (fe_data.update_each & update_gradients &&
131  cell_similarity != CellSimilarity::translation)
132  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
133  mapping.transform(make_array_view(fe_data.shape_gradients, k),
135  mapping_internal,
136  make_array_view(output_data.shape_gradients, k));
137 
138  if (fe_data.update_each & update_hessians &&
139  cell_similarity != CellSimilarity::translation)
140  {
141  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
142  mapping.transform(make_array_view(fe_data.shape_hessians, k),
144  mapping_internal,
145  make_array_view(output_data.shape_hessians, k));
146 
147  if (need_to_correct_higher_derivatives)
148  correct_hessians(output_data, mapping_data, quadrature.size());
149  }
150 
151  if (fe_data.update_each & update_3rd_derivatives &&
152  cell_similarity != CellSimilarity::translation)
153  {
154  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
155  mapping.transform(make_array_view(fe_data.shape_3rd_derivatives, k),
157  mapping_internal,
159  k));
160 
161  if (need_to_correct_higher_derivatives)
162  correct_third_derivatives(output_data, mapping_data, quadrature.size());
163  }
164 }
165 
166 
167 #include "fe_poly.inst"
168 
polynomials_rannacher_turek.h
fe_values.h
tensor_product_polynomials.h
Mapping::transform
virtual void transform(const ArrayView< const Tensor< 1, dim >> &input, const MappingKind kind, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const ArrayView< Tensor< 1, spacedim >> &output) const =0
internal::FEValuesImplementation::FiniteElementRelatedData::shape_gradients
GradientVector shape_gradients
Definition: fe_update_flags.h:590
internal::FEValuesImplementation::FiniteElementRelatedData::shape_hessians
HessianVector shape_hessians
Definition: fe_update_flags.h:597
internal::FEValuesImplementation::FiniteElementRelatedData
Definition: fe_update_flags.h:524
tensor_product_polynomials_bubbles.h
mapping_covariant_gradient
@ mapping_covariant_gradient
Definition: mapping.h:84
update_3rd_derivatives
@ update_3rd_derivatives
Third derivatives of shape functions.
Definition: fe_update_flags.h:96
CellSimilarity::translation
@ translation
Definition: fe_update_flags.h:388
fe_poly.h
CellSimilarity::Similarity
Similarity
Definition: fe_update_flags.h:378
polynomials_piecewise.h
internal::FEValuesImplementation::FiniteElementRelatedData::shape_3rd_derivatives
ThirdDerivativeVector shape_3rd_derivatives
Definition: fe_update_flags.h:604
Mapping
Abstract base class for mapping classes.
Definition: mapping.h:302
mapping_covariant_hessian
@ mapping_covariant_hessian
Definition: mapping.h:134
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:358
update_gradients
@ update_gradients
Shape function gradients.
Definition: fe_update_flags.h:84
ArrayView::make_array_view
ArrayView< typename std::remove_reference< typename std::iterator_traits< Iterator >::reference >::type, MemorySpaceType > make_array_view(const Iterator begin, const Iterator end)
Definition: array_view.h:607
StandardExceptions::ExcInternalError
static ::ExceptionBase & ExcInternalError()
Assert
#define Assert(cond, exc)
Definition: exceptions.h:1419
update_hessians
@ update_hessians
Second derivatives of shape functions.
Definition: fe_update_flags.h:90
Mapping::InternalDataBase
Definition: mapping.h:597
tensor_product_polynomials_const.h
FE_Poly
Definition: fe_poly.h:76
FiniteElement::InternalDataBase
Definition: fe.h:682
Quadrature::size
unsigned int size() const
mapping_covariant
@ mapping_covariant
Definition: mapping.h:73
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:359
Quadrature< 1 >
polynomials_p.h
TriaIterator
Definition: tria_iterator.h:578
qprojector.h