Reference documentation for deal.II version GIT c00406db16 2023-09-28 18:00: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\}}\)
fe_update_flags.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2023 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 
22 #include <deal.II/base/table.h>
23 #include <deal.II/base/tensor.h>
24 
25 #include <vector>
26 
27 
29 
30 // Forward declaration
31 #ifndef DOXYGEN
32 template <int, int>
33 class FiniteElement;
34 #endif
35 
65 {
69 
76  update_values = 0x0001,
78 
82  update_gradients = 0x0002,
84 
88  update_hessians = 0x0004,
90 
96 
103 
130 
137 
144 
150 
155 
161 
168 
175 
181 
187 
217  update_rescale = 0x2000000,
219 
227  // Direct data
235  // Transformation dependence
238  // Volume data
240  // Hermite needs several DOFs to be rescaled each time
242 };
243 
244 
250 template <typename StreamType>
251 inline StreamType &
252 operator<<(StreamType &s, const UpdateFlags u)
253 {
254  s << " UpdateFlags|";
255  if (u & update_values)
256  s << "values|";
257  if (u & update_gradients)
258  s << "gradients|";
259  if (u & update_hessians)
260  s << "hessians|";
261  if (u & update_3rd_derivatives)
262  s << "3rd_derivatives|";
263  if (u & update_quadrature_points)
264  s << "quadrature_points|";
265  if (u & update_JxW_values)
266  s << "JxW_values|";
267  if (u & update_normal_vectors)
268  s << "normal_vectors|";
269  if (u & update_jacobians)
270  s << "jacobians|";
271  if (u & update_inverse_jacobians)
272  s << "inverse_jacobians|";
273  if (u & update_jacobian_grads)
274  s << "jacobian_grads|";
276  s << "covariant_transformation|";
278  s << "contravariant_transformation|";
280  s << "transformation_values|";
282  s << "transformation_gradients|";
284  s << "jacobian_pushed_forward_grads|";
286  s << "jacobian_2nd_derivatives|";
288  s << "jacobian_pushed_forward_2nd_derivatives|";
290  s << "jacobian_3rd_derivatives|";
292  s << "jacobian_pushed_forward_3rd_derivatives|";
293 
294  // TODO: check that 'u' really only has the flags set that are handled above
295  return s;
296 }
297 
298 
308 inline UpdateFlags
309 operator|(const UpdateFlags f1, const UpdateFlags f2)
310 {
311  return static_cast<UpdateFlags>(static_cast<unsigned int>(f1) |
312  static_cast<unsigned int>(f2));
313 }
314 
315 
316 
323 inline UpdateFlags &
325 {
326  f1 = f1 | f2;
327  return f1;
328 }
329 
330 
340 inline UpdateFlags
341 operator&(const UpdateFlags f1, const UpdateFlags f2)
342 {
343  return static_cast<UpdateFlags>(static_cast<unsigned int>(f1) &
344  static_cast<unsigned int>(f2));
345 }
346 
347 
354 inline UpdateFlags &
356 {
357  f1 = f1 & f2;
358  return f1;
359 }
360 
361 
362 
373 namespace CellSimilarity
374 {
376  {
394  };
395 }
396 
397 
398 namespace internal
399 {
400  namespace FEValuesImplementation
401  {
410  template <int dim, int spacedim = dim>
412  {
413  public:
417  void
418  initialize(const unsigned int n_quadrature_points,
420  const UpdateFlags flags);
421 
426  std::size_t
427  memory_consumption() const;
428 
448 
454 
459 
464 
471 
478 
485 
492 
522  std::vector<unsigned int> shape_function_to_row_table;
523  };
524  } // namespace FEValuesImplementation
525 } // namespace internal
526 
527 
533 
534 #endif
OutputOperator< VectorType > & operator<<(OutputOperator< VectorType > &out, unsigned int step)
Definition: operator.h:165
void initialize(const unsigned int n_quadrature_points, const FiniteElement< dim, spacedim > &fe, const UpdateFlags flags)
Definition: fe_values.cc:84
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
IteratorRange< FilteredIterator< BaseIterator > > operator|(IteratorRange< BaseIterator > i, const Predicate &p)
UpdateFlags & operator|=(UpdateFlags &f1, const UpdateFlags f2)
UpdateFlags operator&(const UpdateFlags f1, const UpdateFlags f2)
UpdateFlags
UpdateFlags & operator&=(UpdateFlags &f1, const UpdateFlags f2)
@ update_jacobian_pushed_forward_2nd_derivatives
@ update_volume_elements
Determinant of the Jacobian.
@ update_contravariant_transformation
Contravariant transformation.
@ update_jacobian_pushed_forward_grads
@ update_rescale
@ update_hessians
Second derivatives of shape functions.
@ update_jacobian_3rd_derivatives
@ update_values
Shape function values.
@ update_jacobian_grads
Gradient of volume element.
@ update_normal_vectors
Normal vectors.
@ update_3rd_derivatives
Third derivatives of shape functions.
@ update_JxW_values
Transformed quadrature weights.
@ update_covariant_transformation
Covariant transformation.
@ update_transformation_gradients
Shape function gradients of transformation.
@ update_jacobians
Volume element.
@ update_mapping
@ update_inverse_jacobians
Volume element.
@ update_gradients
Shape function gradients.
@ update_quadrature_points
Transformed quadrature points.
@ update_default
No update.
@ update_jacobian_pushed_forward_3rd_derivatives
@ update_transformation_values
Shape function values of transformation.
@ update_piola
Values needed for Piola transform.
@ update_boundary_forms
Outer normal vector, not normalized.
@ update_jacobian_2nd_derivatives