Reference documentation for deal.II version Git 7026f387cc 2019-10-15 14:19:01 -0400
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
fe_update_flags.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 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 #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/derivative_form.h>
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 
30 DEAL_II_NAMESPACE_OPEN
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 
146 
152 
157 
163 
170 
177 
183 
189 
225 
233  // Direct data
241  // Transformation dependence
244  // Volume data
246 };
247 
248 
254 template <class StreamType>
255 inline StreamType &
256 operator<<(StreamType &s, const UpdateFlags u)
257 {
258  s << " UpdateFlags|";
259  if (u & update_values)
260  s << "values|";
261  if (u & update_gradients)
262  s << "gradients|";
263  if (u & update_hessians)
264  s << "hessians|";
265  if (u & update_3rd_derivatives)
266  s << "3rd_derivatives|";
267  if (u & update_quadrature_points)
268  s << "quadrature_points|";
269  if (u & update_JxW_values)
270  s << "JxW_values|";
271  if (u & update_normal_vectors)
272  s << "normal_vectors|";
273  if (u & update_jacobians)
274  s << "jacobians|";
275  if (u & update_inverse_jacobians)
276  s << "inverse_jacobians|";
277  if (u & update_jacobian_grads)
278  s << "jacobian_grads|";
280  s << "covariant_transformation|";
282  s << "contravariant_transformation|";
284  s << "transformation_values|";
286  s << "transformation_gradients|";
288  s << "jacobian_pushed_forward_grads|";
290  s << "jacobian_2nd_derivatives|";
292  s << "jacobian_pushed_forward_2nd_derivatives|";
294  s << "jacobian_3rd_derivatives|";
296  s << "jacobian_pushed_forward_3rd_derivatives|";
297 
298  // TODO: check that 'u' really only has the flags set that are handled above
299  return s;
300 }
301 
302 
312 inline UpdateFlags
313 operator|(const UpdateFlags f1, const UpdateFlags f2)
314 {
315  return static_cast<UpdateFlags>(static_cast<unsigned int>(f1) |
316  static_cast<unsigned int>(f2));
317 }
318 
319 
320 
327 inline UpdateFlags &
329 {
330  f1 = f1 | f2;
331  return f1;
332 }
333 
334 
344 inline UpdateFlags operator&(const UpdateFlags f1, const UpdateFlags f2)
345 {
346  return static_cast<UpdateFlags>(static_cast<unsigned int>(f1) &
347  static_cast<unsigned int>(f2));
348 }
349 
350 
357 inline UpdateFlags &
359 {
360  f1 = f1 & f2;
361  return f1;
362 }
363 
364 
365 
376 namespace CellSimilarity
377 {
379  {
397  };
398 }
399 
400 
401 namespace internal
402 {
403  namespace FEValuesImplementation
404  {
417  template <int dim, int spacedim = dim>
419  {
420  public:
424  void
425  initialize(const unsigned int n_quadrature_points,
426  const UpdateFlags flags);
427 
432  std::size_t
433  memory_consumption() const;
434 
448  std::vector<double> JxW_values;
449 
453  std::vector<DerivativeForm<1, dim, spacedim>> jacobians;
454 
459  std::vector<DerivativeForm<2, dim, spacedim>> jacobian_grads;
460 
464  std::vector<DerivativeForm<1, spacedim, dim>> inverse_jacobians;
465 
470  std::vector<Tensor<3, spacedim>> jacobian_pushed_forward_grads;
471 
476  std::vector<DerivativeForm<3, dim, spacedim>> jacobian_2nd_derivatives;
477 
482  std::vector<Tensor<4, spacedim>> jacobian_pushed_forward_2nd_derivatives;
483 
488  std::vector<DerivativeForm<4, dim, spacedim>> jacobian_3rd_derivatives;
489 
494  std::vector<Tensor<5, spacedim>> jacobian_pushed_forward_3rd_derivatives;
495 
501  std::vector<Point<spacedim>> quadrature_points;
502 
506  std::vector<Tensor<1, spacedim>> normal_vectors;
507 
511  std::vector<Tensor<1, spacedim>> boundary_forms;
512  };
513 
514 
523  template <int dim, int spacedim = dim>
525  {
526  public:
530  void
531  initialize(const unsigned int n_quadrature_points,
533  const UpdateFlags flags);
534 
539  std::size_t
540  memory_consumption() const;
541 
561 
567 
572 
577 
584 
591 
598 
605 
635  std::vector<unsigned int> shape_function_to_row_table;
636  };
637  } // namespace FEValuesImplementation
638 } // namespace internal
639 
640 
645 DEAL_II_NAMESPACE_CLOSE
646 
647 #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.
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:170
UpdateFlags & operator|=(UpdateFlags &f1, const UpdateFlags f2)
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