Reference documentation for deal.II version Git 5a2787e538 2021-09-21 14:55:10 -0600
\(\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_update_flags.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 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_update_flags_h
17 #define dealii_fe_update_flags_h
18 
19 
20 #include <deal.II/base/config.h>
21 
23 #include <deal.II/base/point.h>
24 #include <deal.II/base/table.h>
25 #include <deal.II/base/tensor.h>
26 
27 #include <vector>
28 
29 
31 
32 // Forward declaration
33 #ifndef DOXYGEN
34 template <int, int>
35 class FiniteElement;
36 #endif
37 
40 
67 {
71 
78  update_values = 0x0001,
80 
84  update_gradients = 0x0002,
86 
90  update_hessians = 0x0004,
92 
98 
105 
124 
131 
138 
144 
149 
155 
162 
169 
175 
181 
209 
217  // Direct data
225  // Transformation dependence
228  // Volume data
230 };
231 
232 
238 template <class StreamType>
239 inline StreamType &
240 operator<<(StreamType &s, const UpdateFlags u)
241 {
242  s << " UpdateFlags|";
243  if (u & update_values)
244  s << "values|";
245  if (u & update_gradients)
246  s << "gradients|";
247  if (u & update_hessians)
248  s << "hessians|";
249  if (u & update_3rd_derivatives)
250  s << "3rd_derivatives|";
251  if (u & update_quadrature_points)
252  s << "quadrature_points|";
253  if (u & update_JxW_values)
254  s << "JxW_values|";
255  if (u & update_normal_vectors)
256  s << "normal_vectors|";
257  if (u & update_jacobians)
258  s << "jacobians|";
259  if (u & update_inverse_jacobians)
260  s << "inverse_jacobians|";
261  if (u & update_jacobian_grads)
262  s << "jacobian_grads|";
264  s << "covariant_transformation|";
266  s << "contravariant_transformation|";
268  s << "transformation_values|";
270  s << "transformation_gradients|";
272  s << "jacobian_pushed_forward_grads|";
274  s << "jacobian_2nd_derivatives|";
276  s << "jacobian_pushed_forward_2nd_derivatives|";
278  s << "jacobian_3rd_derivatives|";
280  s << "jacobian_pushed_forward_3rd_derivatives|";
281 
282  // TODO: check that 'u' really only has the flags set that are handled above
283  return s;
284 }
285 
286 
296 inline UpdateFlags
297 operator|(const UpdateFlags f1, const UpdateFlags f2)
298 {
299  return static_cast<UpdateFlags>(static_cast<unsigned int>(f1) |
300  static_cast<unsigned int>(f2));
301 }
302 
303 
304 
311 inline UpdateFlags &
313 {
314  f1 = f1 | f2;
315  return f1;
316 }
317 
318 
328 inline UpdateFlags
329 operator&(const UpdateFlags f1, const UpdateFlags f2)
330 {
331  return static_cast<UpdateFlags>(static_cast<unsigned int>(f1) &
332  static_cast<unsigned int>(f2));
333 }
334 
335 
342 inline UpdateFlags &
344 {
345  f1 = f1 & f2;
346  return f1;
347 }
348 
349 
350 
361 namespace CellSimilarity
362 {
364  {
382  };
383 }
384 
385 
386 namespace internal
387 {
388  namespace FEValuesImplementation
389  {
402  template <int dim, int spacedim = dim>
404  {
405  public:
409  void
410  initialize(const unsigned int n_quadrature_points,
411  const UpdateFlags flags);
412 
417  std::size_t
418  memory_consumption() const;
419 
433  std::vector<double> JxW_values;
434 
438  std::vector<DerivativeForm<1, dim, spacedim>> jacobians;
439 
444  std::vector<DerivativeForm<2, dim, spacedim>> jacobian_grads;
445 
449  std::vector<DerivativeForm<1, spacedim, dim>> inverse_jacobians;
450 
455  std::vector<Tensor<3, spacedim>> jacobian_pushed_forward_grads;
456 
461  std::vector<DerivativeForm<3, dim, spacedim>> jacobian_2nd_derivatives;
462 
467  std::vector<Tensor<4, spacedim>> jacobian_pushed_forward_2nd_derivatives;
468 
473  std::vector<DerivativeForm<4, dim, spacedim>> jacobian_3rd_derivatives;
474 
479  std::vector<Tensor<5, spacedim>> jacobian_pushed_forward_3rd_derivatives;
480 
486  std::vector<Point<spacedim>> quadrature_points;
487 
491  std::vector<Tensor<1, spacedim>> normal_vectors;
492 
496  std::vector<Tensor<1, spacedim>> boundary_forms;
497  };
498 
499 
508  template <int dim, int spacedim = dim>
510  {
511  public:
515  void
516  initialize(const unsigned int n_quadrature_points,
518  const UpdateFlags flags);
519 
524  std::size_t
525  memory_consumption() const;
526 
546 
552 
557 
562 
569 
576 
583 
590 
620  std::vector<unsigned int> shape_function_to_row_table;
621  };
622  } // namespace FEValuesImplementation
623 } // namespace internal
624 
625 
631 
632 #endif
Transformed quadrature weights.
Shape function values.
Contravariant transformation.
std::vector< Tensor< 1, spacedim > > boundary_forms
Volume element.
std::vector< DerivativeForm< 3, dim, spacedim > > jacobian_2nd_derivatives
Outer normal vector, not normalized.
Determinant of the Jacobian.
std::vector< DerivativeForm< 4, dim, spacedim > > jacobian_3rd_derivatives
Transformed quadrature points.
std::vector< Tensor< 4, spacedim > > jacobian_pushed_forward_2nd_derivatives
UpdateFlags operator|(const UpdateFlags f1, const UpdateFlags f2)
Shape function gradients of transformation.
UpdateFlags operator &(const UpdateFlags f1, const UpdateFlags f2)
No update.
Third derivatives of shape functions.
UpdateFlags
Shape function values of transformation.
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
UpdateFlags & operator &=(UpdateFlags &f1, const UpdateFlags f2)
Second derivatives of shape functions.
Gradient of volume element.
std::vector< Tensor< 3, spacedim > > jacobian_pushed_forward_grads
std::vector< Point< spacedim > > quadrature_points
OutputOperator< VectorType > & operator<<(OutputOperator< VectorType > &out, unsigned int step)
Definition: operator.h:165
UpdateFlags & operator|=(UpdateFlags &f1, const UpdateFlags f2)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
Shape function gradients.
Normal vectors.
std::vector< DerivativeForm< 1, dim, spacedim > > jacobians
std::vector< Tensor< 5, spacedim > > jacobian_pushed_forward_3rd_derivatives
std::vector< DerivativeForm< 2, dim, spacedim > > jacobian_grads
std::vector< DerivativeForm< 1, spacedim, dim > > inverse_jacobians
Values needed for Piola transform.
Covariant transformation.
std::vector< Tensor< 1, spacedim > > normal_vectors
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)