Loading [MathJax]/extensions/TeX/newcommand.js
 Reference documentation for deal.II version 9.1.1
\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
fe_poly.cc
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2004 - 2018 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 
17 #include <deal.II/base/polynomials_p.h>
18 #include <deal.II/base/polynomials_piecewise.h>
19 #include <deal.II/base/polynomials_rannacher_turek.h>
20 #include <deal.II/base/qprojector.h>
21 #include <deal.II/base/tensor_product_polynomials.h>
22 #include <deal.II/base/tensor_product_polynomials_bubbles.h>
23 #include <deal.II/base/tensor_product_polynomials_const.h>
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 
29 DEAL_II_NAMESPACE_OPEN
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  // transform gradients and higher derivatives. there is nothing to do
54  // for values since we already emplaced them into output_data when
55  // we were in get_data()
56  if (fe_data.update_each & update_gradients &&
57  cell_similarity != CellSimilarity::translation)
58  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
59  mapping.transform(make_array_view(fe_data.shape_gradients, k),
61  mapping_internal,
62  make_array_view(output_data.shape_gradients, k));
63 
64  if (fe_data.update_each & update_hessians &&
65  cell_similarity != CellSimilarity::translation)
66  {
67  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
68  mapping.transform(make_array_view(fe_data.shape_hessians, k),
70  mapping_internal,
71  make_array_view(output_data.shape_hessians, k));
72 
73  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
74  for (unsigned int i = 0; i < quadrature.size(); ++i)
75  for (unsigned int j = 0; j < 2; ++j)
76  output_data.shape_hessians[k][i] -=
77  mapping_data.jacobian_pushed_forward_grads[i][j] *
78  output_data.shape_gradients[k][i][j];
79  }
80 
81  if (fe_data.update_each & update_3rd_derivatives &&
82  cell_similarity != CellSimilarity::translation)
83  {
84  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
85  mapping.transform(make_array_view(fe_data.shape_3rd_derivatives, k),
87  mapping_internal,
89  k));
90 
91  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
92  correct_third_derivatives(output_data,
93  mapping_data,
94  quadrature.size(),
95  k);
96  }
97 }
98 
99 
100 
101 template <>
102 void
103 FE_Poly<TensorProductPolynomials<2>, 2, 3>::fill_fe_values(
105  const CellSimilarity::Similarity cell_similarity,
106  const Quadrature<2> & quadrature,
107  const Mapping<2, 3> & mapping,
108  const Mapping<2, 3>::InternalDataBase &mapping_internal,
109  const ::internal::FEValuesImplementation::MappingRelatedData<2, 3>
110  & mapping_data,
111  const FiniteElement<2, 3>::InternalDataBase &fe_internal,
113  &output_data) const
114 {
115  // convert data object to internal data for this class. fails with an
116  // exception if that is not possible
117  Assert(dynamic_cast<const InternalData *>(&fe_internal) != nullptr,
118  ExcInternalError());
119  const InternalData &fe_data =
120  static_cast<const InternalData &>(fe_internal); // NOLINT
121 
122  // transform gradients and higher derivatives. there is nothing to do
123  // for values since we already emplaced them into output_data when
124  // we were in get_data()
125  if (fe_data.update_each & update_gradients &&
126  cell_similarity != CellSimilarity::translation)
127  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
128  mapping.transform(make_array_view(fe_data.shape_gradients, k),
130  mapping_internal,
131  make_array_view(output_data.shape_gradients, k));
132 
133  if (fe_data.update_each & update_hessians &&
134  cell_similarity != CellSimilarity::translation)
135  {
136  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
137  mapping.transform(make_array_view(fe_data.shape_hessians, k),
139  mapping_internal,
140  make_array_view(output_data.shape_hessians, k));
141 
142  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
143  for (unsigned int i = 0; i < quadrature.size(); ++i)
144  for (unsigned int j = 0; j < 3; ++j)
145  output_data.shape_hessians[k][i] -=
146  mapping_data.jacobian_pushed_forward_grads[i][j] *
147  output_data.shape_gradients[k][i][j];
148  }
149 
150  if (fe_data.update_each & update_3rd_derivatives &&
151  cell_similarity != CellSimilarity::translation)
152  {
153  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
154  mapping.transform(make_array_view(fe_data.shape_3rd_derivatives, k),
156  mapping_internal,
158  k));
159 
160  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
161  correct_third_derivatives(output_data,
162  mapping_data,
163  quadrature.size(),
164  k);
165  }
166 }
167 
168 
169 template <>
170 void
171 FE_Poly<PolynomialSpace<1>, 1, 2>::fill_fe_values(
173  const CellSimilarity::Similarity cell_similarity,
174  const Quadrature<1> & quadrature,
175  const Mapping<1, 2> & mapping,
176  const Mapping<1, 2>::InternalDataBase &mapping_internal,
177  const ::internal::FEValuesImplementation::MappingRelatedData<1, 2>
178  & mapping_data,
179  const FiniteElement<1, 2>::InternalDataBase &fe_internal,
181  &output_data) const
182 {
183  // convert data object to internal data for this class. fails with an
184  // exception if that is not possible
185  Assert(dynamic_cast<const InternalData *>(&fe_internal) != nullptr,
186  ExcInternalError());
187  const InternalData &fe_data =
188  static_cast<const InternalData &>(fe_internal); // NOLINT
189 
190  // transform gradients and higher derivatives. there is nothing to do
191  // for values since we already emplaced them into output_data when
192  // we were in get_data()
193  if (fe_data.update_each & update_gradients &&
194  cell_similarity != CellSimilarity::translation)
195  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
196  mapping.transform(make_array_view(fe_data.shape_gradients, k),
198  mapping_internal,
199  make_array_view(output_data.shape_gradients, k));
200 
201  if (fe_data.update_each & update_hessians &&
202  cell_similarity != CellSimilarity::translation)
203  {
204  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
205  mapping.transform(make_array_view(fe_data.shape_hessians, k),
207  mapping_internal,
208  make_array_view(output_data.shape_hessians, k));
209 
210  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
211  for (unsigned int i = 0; i < quadrature.size(); ++i)
212  for (unsigned int j = 0; j < 2; ++j)
213  output_data.shape_hessians[k][i] -=
214  mapping_data.jacobian_pushed_forward_grads[i][j] *
215  output_data.shape_gradients[k][i][j];
216  }
217 
218  if (fe_data.update_each & update_3rd_derivatives &&
219  cell_similarity != CellSimilarity::translation)
220  {
221  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
222  mapping.transform(make_array_view(fe_data.shape_3rd_derivatives, k),
224  mapping_internal,
226  k));
227 
228  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
229  correct_third_derivatives(output_data,
230  mapping_data,
231  quadrature.size(),
232  k);
233  }
234 }
235 
236 
237 template <>
238 void
239 FE_Poly<PolynomialSpace<2>, 2, 3>::fill_fe_values(
241  const CellSimilarity::Similarity cell_similarity,
242  const Quadrature<2> & quadrature,
243  const Mapping<2, 3> & mapping,
244  const Mapping<2, 3>::InternalDataBase &mapping_internal,
245  const ::internal::FEValuesImplementation::MappingRelatedData<2, 3>
246  & mapping_data,
247  const FiniteElement<2, 3>::InternalDataBase &fe_internal,
249  &output_data) const
250 {
251  // convert data object to internal data for this class. fails with an
252  // exception if that is not possible
253  Assert(dynamic_cast<const InternalData *>(&fe_internal) != nullptr,
254  ExcInternalError());
255  const InternalData &fe_data =
256  static_cast<const InternalData &>(fe_internal); // NOLINT
257 
258  // transform gradients and higher derivatives. there is nothing to do
259  // for values since we already emplaced them into output_data when
260  // we were in get_data()
261  if (fe_data.update_each & update_gradients &&
262  cell_similarity != CellSimilarity::translation)
263  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
264  mapping.transform(make_array_view(fe_data.shape_gradients, k),
266  mapping_internal,
267  make_array_view(output_data.shape_gradients, k));
268 
269  if (fe_data.update_each & update_hessians &&
270  cell_similarity != CellSimilarity::translation)
271  {
272  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
273  mapping.transform(make_array_view(fe_data.shape_hessians, k),
275  mapping_internal,
276  make_array_view(output_data.shape_hessians, k));
277 
278  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
279  for (unsigned int i = 0; i < quadrature.size(); ++i)
280  for (unsigned int j = 0; j < 3; ++j)
281  output_data.shape_hessians[k][i] -=
282  mapping_data.jacobian_pushed_forward_grads[i][j] *
283  output_data.shape_gradients[k][i][j];
284  }
285 
286  if (fe_data.update_each & update_3rd_derivatives &&
287  cell_similarity != CellSimilarity::translation)
288  {
289  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
290  mapping.transform(make_array_view(fe_data.shape_3rd_derivatives, k),
292  mapping_internal,
294  k));
295 
296  for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
297  correct_third_derivatives(output_data,
298  mapping_data,
299  quadrature.size(),
300  k);
301  }
302 }
303 
304 
305 #include "fe_poly.inst"
306 
307 DEAL_II_NAMESPACE_CLOSE
virtual void transform(const ArrayView< const Tensor< 1, dim >> &input, const MappingType type, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const ArrayView< Tensor< 1, spacedim >> &output) const =0
Third derivatives of shape functions.
#define Assert(cond, exc)
Definition: exceptions.h:1407
Abstract base class for mapping classes.
Definition: dof_tools.h:57
Second derivatives of shape functions.
unsigned int size() const
Shape function gradients.
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:594
static ::ExceptionBase & ExcInternalError()