Reference documentation for deal.II version 9.0.0
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 at
12 // the top level of the deal.II distribution.
13 //
14 // ---------------------------------------------------------------------
15 
16 
17 #include <deal.II/base/qprojector.h>
18 #include <deal.II/base/tensor_product_polynomials.h>
19 #include <deal.II/base/tensor_product_polynomials_const.h>
20 #include <deal.II/base/tensor_product_polynomials_bubbles.h>
21 #include <deal.II/base/polynomials_p.h>
22 #include <deal.II/base/polynomials_piecewise.h>
23 #include <deal.II/base/polynomials_rannacher_turek.h>
24 #include <deal.II/fe/fe_poly.h>
25 #include <deal.II/fe/fe_values.h>
26 #include <deal.II/fe/fe_poly.templates.h>
27 
28 DEAL_II_NAMESPACE_OPEN
29 
30 
31 template <>
32 void
34 fill_fe_values (const Triangulation<1,2>::cell_iterator &,
35  const CellSimilarity::Similarity cell_similarity,
36  const Quadrature<1> &quadrature,
37  const Mapping<1,2> &mapping,
38  const Mapping<1,2>::InternalDataBase &mapping_internal,
39  const ::internal::FEValuesImplementation::MappingRelatedData<1,2> &mapping_data,
40  const FiniteElement<1,2>::InternalDataBase &fe_internal,
42 {
43  // convert data object to internal
44  // data for this class. fails with
45  // an exception if that is not
46  // possible
47  Assert (dynamic_cast<const InternalData *> (&fe_internal) != nullptr, ExcInternalError());
48  const InternalData &fe_data = static_cast<const InternalData &> (fe_internal);
49 
50  // transform gradients and higher derivatives. there is nothing to do
51  // for values since we already emplaced them into output_data when
52  // we were in get_data()
53  if (fe_data.update_each & update_gradients && cell_similarity != CellSimilarity::translation)
54  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
55  mapping.transform (make_array_view(fe_data.shape_gradients, k),
57  mapping_internal,
58  make_array_view(output_data.shape_gradients, k));
59 
60  if (fe_data.update_each & update_hessians && cell_similarity != CellSimilarity::translation)
61  {
62  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
63  mapping.transform (make_array_view(fe_data.shape_hessians, k),
65  mapping_internal,
66  make_array_view(output_data.shape_hessians, k));
67 
68  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
69  for (unsigned int i=0; i<quadrature.size(); ++i)
70  for (unsigned int j=0; j<2; ++j)
71  output_data.shape_hessians[k][i] -=
72  mapping_data.jacobian_pushed_forward_grads[i][j]
73  * output_data.shape_gradients[k][i][j];
74  }
75 
76  if (fe_data.update_each & update_3rd_derivatives && cell_similarity != CellSimilarity::translation)
77  {
78  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
79  mapping.transform (make_array_view(fe_data.shape_3rd_derivatives, k),
81  mapping_internal,
82  make_array_view(output_data.shape_3rd_derivatives, k));
83 
84  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
85  correct_third_derivatives(output_data, mapping_data, quadrature.size(), k);
86  }
87 }
88 
89 
90 
91 template <>
92 void
94 fill_fe_values (const Triangulation<2,3>::cell_iterator &,
95  const CellSimilarity::Similarity cell_similarity,
96  const Quadrature<2> &quadrature,
97  const Mapping<2,3> &mapping,
98  const Mapping<2,3>::InternalDataBase &mapping_internal,
99  const ::internal::FEValuesImplementation::MappingRelatedData<2,3> &mapping_data,
100  const FiniteElement<2,3>::InternalDataBase &fe_internal,
102 {
103 
104  // assert that the following dynamics
105  // cast is really well-defined.
106  Assert (dynamic_cast<const InternalData *> (&fe_internal) != nullptr, ExcInternalError());
107  const InternalData &fe_data = static_cast<const InternalData &> (fe_internal);
108 
109  // transform gradients and higher derivatives. there is nothing to do
110  // for values since we already emplaced them into output_data when
111  // we were in get_data()
112  if (fe_data.update_each & update_gradients && cell_similarity != CellSimilarity::translation)
113  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
114  mapping.transform (make_array_view(fe_data.shape_gradients, k),
116  mapping_internal,
117  make_array_view(output_data.shape_gradients, k));
118 
119  if (fe_data.update_each & update_hessians && cell_similarity != CellSimilarity::translation)
120  {
121  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
122  mapping.transform (make_array_view(fe_data.shape_hessians, k),
124  mapping_internal,
125  make_array_view(output_data.shape_hessians, k));
126 
127  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
128  for (unsigned int i=0; i<quadrature.size(); ++i)
129  for (unsigned int j=0; j<3; ++j)
130  output_data.shape_hessians[k][i] -=
131  mapping_data.jacobian_pushed_forward_grads[i][j]
132  * output_data.shape_gradients[k][i][j];
133  }
134 
135  if (fe_data.update_each & update_3rd_derivatives && cell_similarity != CellSimilarity::translation)
136  {
137  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
138  mapping.transform (make_array_view(fe_data.shape_3rd_derivatives, k),
140  mapping_internal,
141  make_array_view(output_data.shape_3rd_derivatives, k));
142 
143  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
144  correct_third_derivatives(output_data, mapping_data, quadrature.size(), k);
145  }
146 }
147 
148 
149 template <>
150 void
152 fill_fe_values (const Triangulation<1,2>::cell_iterator &,
153  const CellSimilarity::Similarity cell_similarity,
154  const Quadrature<1> &quadrature,
155  const Mapping<1,2> &mapping,
156  const Mapping<1,2>::InternalDataBase &mapping_internal,
157  const ::internal::FEValuesImplementation::MappingRelatedData<1,2> &mapping_data,
158  const FiniteElement<1,2>::InternalDataBase &fe_internal,
160 {
161  // convert data object to internal
162  // data for this class. fails with
163  // an exception if that is not
164  // possible
165 
166  Assert (dynamic_cast<const InternalData *> (&fe_internal) != nullptr, ExcInternalError());
167  const InternalData &fe_data = static_cast<const InternalData &> (fe_internal);
168 
169  // transform gradients and higher derivatives. there is nothing to do
170  // for values since we already emplaced them into output_data when
171  // we were in get_data()
172  if (fe_data.update_each & update_gradients && cell_similarity != CellSimilarity::translation)
173  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
174  mapping.transform (make_array_view(fe_data.shape_gradients, k),
176  mapping_internal,
177  make_array_view(output_data.shape_gradients, k));
178 
179  if (fe_data.update_each & update_hessians && cell_similarity != CellSimilarity::translation)
180  {
181  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
182  mapping.transform (make_array_view(fe_data.shape_hessians, k),
184  mapping_internal,
185  make_array_view(output_data.shape_hessians, k));
186 
187  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
188  for (unsigned int i=0; i<quadrature.size(); ++i)
189  for (unsigned int j=0; j<2; ++j)
190  output_data.shape_hessians[k][i] -=
191  mapping_data.jacobian_pushed_forward_grads[i][j]
192  * output_data.shape_gradients[k][i][j];
193  }
194 
195  if (fe_data.update_each & update_3rd_derivatives && cell_similarity != CellSimilarity::translation)
196  {
197  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
198  mapping.transform (make_array_view(fe_data.shape_3rd_derivatives, k),
200  mapping_internal,
201  make_array_view(output_data.shape_3rd_derivatives, k));
202 
203  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
204  correct_third_derivatives(output_data, mapping_data, quadrature.size(), k);
205  }
206 }
207 
208 
209 template <>
210 void
212 fill_fe_values (const Triangulation<2,3>::cell_iterator &,
213  const CellSimilarity::Similarity cell_similarity,
214  const Quadrature<2> &quadrature,
215  const Mapping<2,3> &mapping,
216  const Mapping<2,3>::InternalDataBase &mapping_internal,
217  const ::internal::FEValuesImplementation::MappingRelatedData<2,3> &mapping_data,
218  const FiniteElement<2,3>::InternalDataBase &fe_internal,
220 {
221  Assert (dynamic_cast<const InternalData *> (&fe_internal) != nullptr, ExcInternalError());
222  const InternalData &fe_data = static_cast<const InternalData &> (fe_internal);
223 
224  // transform gradients and higher derivatives. there is nothing to do
225  // for values since we already emplaced them into output_data when
226  // we were in get_data()
227  if (fe_data.update_each & update_gradients && cell_similarity != CellSimilarity::translation)
228  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
229  mapping.transform (make_array_view(fe_data.shape_gradients, k),
231  mapping_internal,
232  make_array_view(output_data.shape_gradients, k));
233 
234  if (fe_data.update_each & update_hessians && cell_similarity != CellSimilarity::translation)
235  {
236  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
237  mapping.transform (make_array_view(fe_data.shape_hessians, k),
239  mapping_internal,
240  make_array_view(output_data.shape_hessians, k));
241 
242  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
243  for (unsigned int i=0; i<quadrature.size(); ++i)
244  for (unsigned int j=0; j<3; ++j)
245  output_data.shape_hessians[k][i] -=
246  mapping_data.jacobian_pushed_forward_grads[i][j]
247  * output_data.shape_gradients[k][i][j];
248  }
249 
250  if (fe_data.update_each & update_3rd_derivatives && cell_similarity != CellSimilarity::translation)
251  {
252  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
253  mapping.transform (make_array_view(fe_data.shape_3rd_derivatives, k),
255  mapping_internal,
256  make_array_view(output_data.shape_3rd_derivatives, k));
257 
258  for (unsigned int k=0; k<this->dofs_per_cell; ++k)
259  correct_third_derivatives(output_data, mapping_data, quadrature.size(), k);
260  }
261 }
262 
263 
264 #include "fe_poly.inst"
265 
266 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
ArrayView< typename std::remove_reference< typename std::iterator_traits< Iterator >::reference >::type > make_array_view(const Iterator begin, const Iterator end)
Definition: array_view.h:490
Third derivatives of shape functions.
#define Assert(cond, exc)
Definition: exceptions.h:1142
Abstract base class for mapping classes.
Definition: dof_tools.h:46
Second derivatives of shape functions.
unsigned int size() const
Shape function gradients.
static ::ExceptionBase & ExcInternalError()