Reference documentation for deal.II version GIT relicensing-437-g81ec864850 2024-04-19 07:30:02+00:00
\(\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\}}\)
Loading...
Searching...
No Matches
fe_poly_tensor.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2005 - 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_fe_poly_tensor_h
16#define dealii_fe_poly_tensor_h
17
18
19#include <deal.II/base/config.h>
20
22#include <deal.II/base/mutex.h>
25
26#include <deal.II/fe/fe.h>
27
29
30#include <memory>
31
33
139template <int dim, int spacedim = dim>
140class FE_PolyTensor : public FiniteElement<dim, spacedim>
141{
142public:
150 const FiniteElementData<dim> &fe_data,
151 const std::vector<bool> &restriction_is_additive_flags,
152 const std::vector<ComponentMask> &nonzero_components);
153
154
159
160 // for documentation, see the FiniteElement base class
161 virtual UpdateFlags
162 requires_update_flags(const UpdateFlags update_flags) const override;
163
170 virtual double
171 shape_value(const unsigned int i, const Point<dim> &p) const override;
172
173 // documentation inherited from the base class
174 virtual double
175 shape_value_component(const unsigned int i,
176 const Point<dim> &p,
177 const unsigned int component) const override;
178
185 virtual Tensor<1, dim>
186 shape_grad(const unsigned int i, const Point<dim> &p) const override;
187
188 // documentation inherited from the base class
189 virtual Tensor<1, dim>
190 shape_grad_component(const unsigned int i,
191 const Point<dim> &p,
192 const unsigned int component) const override;
193
200 virtual Tensor<2, dim>
201 shape_grad_grad(const unsigned int i, const Point<dim> &p) const override;
202
203 // documentation inherited from the base class
204 virtual Tensor<2, dim>
205 shape_grad_grad_component(const unsigned int i,
206 const Point<dim> &p,
207 const unsigned int component) const override;
208
209protected:
217 std::vector<MappingKind> mapping_kind;
218
223 bool
225
239 bool
241 const unsigned int face_no,
242 const bool face_orientation,
243 const bool face_flip,
244 const bool face_rotation) const;
245
265
270 get_mapping_kind(const unsigned int i) const;
271
272 /* NOTE: The following function has its definition inlined into the class
273 declaration because we otherwise run into a compiler error with MS Visual
274 Studio. */
275 virtual std::unique_ptr<
278 const UpdateFlags update_flags,
279 const Mapping<dim, spacedim> &mapping,
280 const Quadrature<dim> &quadrature,
282 spacedim>
283 &output_data) const override
284 {
285 (void)mapping;
286 (void)output_data;
287 // generate a new data object and
288 // initialize some fields
289 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
290 data_ptr = std::make_unique<InternalData>();
291 auto &data = dynamic_cast<InternalData &>(*data_ptr);
292 data.update_each = requires_update_flags(update_flags);
293
294 const unsigned int n_q_points = quadrature.size();
295
296 // some scratch arrays
297 std::vector<Tensor<1, dim>> values(0);
298 std::vector<Tensor<2, dim>> grads(0);
299 std::vector<Tensor<3, dim>> grad_grads(0);
300 std::vector<Tensor<4, dim>> third_derivatives(0);
301 std::vector<Tensor<5, dim>> fourth_derivatives(0);
302
303 if (update_flags & (update_values | update_gradients | update_hessians))
304 data.dof_sign_change.resize(this->dofs_per_cell);
305
306 // initialize fields only if really
307 // necessary. otherwise, don't
308 // allocate memory
309
310 const bool update_transformed_shape_values =
311 std::any_of(this->mapping_kind.begin(),
312 this->mapping_kind.end(),
313 [](const MappingKind t) { return t != mapping_none; });
314
315 const bool update_transformed_shape_grads =
316 std::any_of(this->mapping_kind.begin(),
317 this->mapping_kind.end(),
318 [](const MappingKind t) {
319 return (t == mapping_raviart_thomas || t == mapping_piola ||
320 t == mapping_nedelec || t == mapping_contravariant);
321 });
322
323 const bool update_transformed_shape_hessian_tensors =
324 update_transformed_shape_values;
325
326 if (update_flags & update_values)
327 {
328 values.resize(this->n_dofs_per_cell());
329 data.shape_values.reinit(this->n_dofs_per_cell(), n_q_points);
330 if (update_transformed_shape_values)
331 data.transformed_shape_values.resize(n_q_points);
332 }
333
334 if (update_flags & update_gradients)
335 {
336 grads.resize(this->n_dofs_per_cell());
337 data.shape_grads.reinit(this->n_dofs_per_cell(), n_q_points);
338 data.transformed_shape_grads.resize(n_q_points);
339
340 if (update_transformed_shape_grads)
341 data.untransformed_shape_grads.resize(n_q_points);
342 }
343
344 if (update_flags & update_hessians)
345 {
346 grad_grads.resize(this->n_dofs_per_cell());
347 data.shape_grad_grads.reinit(this->n_dofs_per_cell(), n_q_points);
348 data.transformed_shape_hessians.resize(n_q_points);
349 if (update_transformed_shape_hessian_tensors)
350 data.untransformed_shape_hessian_tensors.resize(n_q_points);
351 }
352
353 // Compute shape function values
354 // and derivatives and hessians on
355 // the reference cell.
356 // Make sure, that for the
357 // node values N_i holds
358 // N_i(v_j)=\delta_ij for all basis
359 // functions v_j
360 if (update_flags & (update_values | update_gradients))
361 for (unsigned int k = 0; k < n_q_points; ++k)
362 {
363 poly_space->evaluate(quadrature.point(k),
364 values,
365 grads,
366 grad_grads,
367 third_derivatives,
368 fourth_derivatives);
369
370 if (update_flags & update_values)
371 {
372 if (inverse_node_matrix.n_cols() == 0)
373 for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
374 data.shape_values[i][k] = values[i];
375 else
376 for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
377 {
378 Tensor<1, dim> add_values;
379 for (unsigned int j = 0; j < this->n_dofs_per_cell(); ++j)
380 add_values += inverse_node_matrix(j, i) * values[j];
381 data.shape_values[i][k] = add_values;
382 }
383 }
384
385 if (update_flags & update_gradients)
386 {
387 if (inverse_node_matrix.n_cols() == 0)
388 for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
389 data.shape_grads[i][k] = grads[i];
390 else
391 for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
392 {
393 Tensor<2, dim> add_grads;
394 for (unsigned int j = 0; j < this->n_dofs_per_cell(); ++j)
395 add_grads += inverse_node_matrix(j, i) * grads[j];
396 data.shape_grads[i][k] = add_grads;
397 }
398 }
399
400 if (update_flags & update_hessians)
401 {
402 if (inverse_node_matrix.n_cols() == 0)
403 for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
404 data.shape_grad_grads[i][k] = grad_grads[i];
405 else
406 for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
407 {
408 Tensor<3, dim> add_grad_grads;
409 for (unsigned int j = 0; j < this->n_dofs_per_cell(); ++j)
410 add_grad_grads +=
411 inverse_node_matrix(j, i) * grad_grads[j];
412 data.shape_grad_grads[i][k] = add_grad_grads;
413 }
414 }
415 }
416 return data_ptr;
417 }
418
419 virtual void
422 const CellSimilarity::Similarity cell_similarity,
423 const Quadrature<dim> &quadrature,
424 const Mapping<dim, spacedim> &mapping,
425 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
427 &mapping_data,
428 const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
430 spacedim>
431 &output_data) const override;
432
433 using FiniteElement<dim, spacedim>::fill_fe_face_values;
434
435 virtual void
438 const unsigned int face_no,
439 const hp::QCollection<dim - 1> &quadrature,
440 const Mapping<dim, spacedim> &mapping,
441 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
443 &mapping_data,
444 const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
446 spacedim>
447 &output_data) const override;
448
449 virtual void
452 const unsigned int face_no,
453 const unsigned int sub_no,
454 const Quadrature<dim - 1> &quadrature,
455 const Mapping<dim, spacedim> &mapping,
456 const typename Mapping<dim, spacedim>::InternalDataBase &mapping_internal,
458 &mapping_data,
459 const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
461 spacedim>
462 &output_data) const override;
463
473 class InternalData : public FiniteElement<dim, spacedim>::InternalDataBase
474 {
475 public:
481
488
495
499 mutable std::vector<double> dof_sign_change;
500 mutable std::vector<Tensor<1, spacedim>> transformed_shape_values;
501 // for shape_gradient computations
502 mutable std::vector<Tensor<2, spacedim>> transformed_shape_grads;
503 mutable std::vector<Tensor<2, dim>> untransformed_shape_grads;
504 // for shape_hessian computations
505 mutable std::vector<Tensor<3, spacedim>> transformed_shape_hessians;
506 mutable std::vector<Tensor<3, dim>> untransformed_shape_hessian_tensors;
507 };
508
509
510
515 const std::unique_ptr<const TensorPolynomialsBase<dim>> poly_space;
516
531
536
543
547 mutable std::vector<Tensor<1, dim>> cached_values;
548
552 mutable std::vector<Tensor<2, dim>> cached_grads;
553
558 mutable std::vector<Tensor<3, dim>> cached_grad_grads;
559};
560
562
563#endif
std::vector< Tensor< 3, spacedim > > transformed_shape_hessians
Table< 2, DerivativeForm< 2, dim, spacedim > > shape_grad_grads
Table< 2, DerivativeForm< 1, dim, spacedim > > shape_grads
std::vector< Tensor< 2, spacedim > > transformed_shape_grads
std::vector< Tensor< 3, dim > > untransformed_shape_hessian_tensors
std::vector< Tensor< 2, dim > > untransformed_shape_grads
Table< 2, Tensor< 1, dim > > shape_values
std::vector< Tensor< 1, spacedim > > transformed_shape_values
std::vector< double > dof_sign_change
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
std::vector< Table< 2, bool > > adjust_quad_dof_sign_for_face_orientation_table
std::vector< Tensor< 2, dim > > cached_grads
FullMatrix< double > inverse_node_matrix
virtual Tensor< 1, dim > shape_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
FE_PolyTensor(const FE_PolyTensor &fe)
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
Point< dim > cached_point
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const override
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const hp::QCollection< 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
Threads::Mutex cache_mutex
virtual std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
virtual double shape_value(const unsigned int i, const Point< dim > &p) const override
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
bool adjust_quad_dof_sign_for_face_orientation(const unsigned int index, const unsigned int face_no, const bool face_orientation, const bool face_flip, const bool face_rotation) const
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
std::vector< MappingKind > mapping_kind
MappingKind get_mapping_kind(const unsigned int i) const
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)
std::vector< Tensor< 1, dim > > cached_values
bool single_mapping_kind() const
const std::unique_ptr< const TensorPolynomialsBase< dim > > poly_space
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
std::vector< Tensor< 3, dim > > cached_grad_grads
unsigned int n_dofs_per_cell() const
const unsigned int dofs_per_cell
Definition fe_data.h:436
const std::vector< bool > restriction_is_additive_flags
Definition fe.h:2564
const std::vector< ComponentMask > nonzero_components
Definition fe.h:2573
Abstract base class for mapping classes.
Definition mapping.h:318
Definition point.h:111
const Point< dim > & point(const unsigned int i) const
unsigned int size() const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
UpdateFlags
@ update_hessians
Second derivatives of shape functions.
@ update_values
Shape function values.
@ update_gradients
Shape function gradients.
MappingKind
Definition mapping.h:79