Reference documentation for deal.II version Git 2618e0f 2017-11-23 17:25:26 +0100
fe_update_flags.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2016 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 #ifndef dealii_fe_update_flags_h
17 #define dealii_fe_update_flags_h
18 
19 
20 #include <deal.II/base/config.h>
21 #include <deal.II/base/table.h>
22 #include <deal.II/base/derivative_form.h>
23 #include <deal.II/base/point.h>
24 #include <deal.II/base/tensor.h>
25 
26 #include <vector>
27 
28 
29 DEAL_II_NAMESPACE_OPEN
30 
31 template <int,int> class FiniteElement;
32 
33 
36 
63 {
67 
74  update_values = 0x0001,
76 
80  update_gradients = 0x0002,
82 
86  update_hessians = 0x0004,
88 
94 
101 
106 
113 
128 
134 
139 
145 
152 
159 
165 
171 
207 
215  // Direct data
228  // Transformation dependence
233  // Volume data
235 };
236 
237 
243 template <class StreamType>
244 inline
245 StreamType &operator << (StreamType &s,
246  const UpdateFlags u)
247 {
248  s << " UpdateFlags|";
249  if (u & update_values) s << "values|";
250  if (u & update_gradients) s << "gradients|";
251  if (u & update_hessians) s << "hessians|";
252  if (u & update_3rd_derivatives) s << "3rd_derivatives|";
253  if (u & update_quadrature_points) s << "quadrature_points|";
254  if (u & update_JxW_values) s << "JxW_values|";
255  if (u & update_normal_vectors) s << "normal_vectors|";
256  if (u & update_jacobians) s << "jacobians|";
257  if (u & update_inverse_jacobians) s << "inverse_jacobians|";
258  if (u & update_jacobian_grads) s << "jacobian_grads|";
259  if (u & update_covariant_transformation) s << "covariant_transformation|";
260  if (u & update_contravariant_transformation) s << "contravariant_transformation|";
261  if (u & update_transformation_values) s << "transformation_values|";
262  if (u & update_transformation_gradients) s << "transformation_gradients|";
263  if (u & update_jacobian_pushed_forward_grads) s << "jacobian_pushed_forward_grads|";
264  if (u & update_jacobian_2nd_derivatives) s << "jacobian_2nd_derivatives|";
265  if (u & update_jacobian_pushed_forward_2nd_derivatives) s << "jacobian_pushed_forward_2nd_derivatives|";
266  if (u &update_jacobian_3rd_derivatives) s << "jacobian_3rd_derivatives|";
267  if (u & update_jacobian_pushed_forward_3rd_derivatives) s << "jacobian_pushed_forward_3rd_derivatives|";
268 
269 //TODO: check that 'u' really only has the flags set that are handled above
270  return s;
271 }
272 
273 
283 inline
286  const UpdateFlags f2)
287 {
288  return static_cast<UpdateFlags> (
289  static_cast<unsigned int> (f1) |
290  static_cast<unsigned int> (f2));
291 }
292 
293 
294 
295 
302 inline
303 UpdateFlags &
305  const UpdateFlags f2)
306 {
307  f1 = f1 | f2;
308  return f1;
309 }
310 
311 
321 inline
324  const UpdateFlags f2)
325 {
326  return static_cast<UpdateFlags> (
327  static_cast<unsigned int> (f1) &
328  static_cast<unsigned int> (f2));
329 }
330 
331 
338 inline
339 UpdateFlags &
341  const UpdateFlags f2)
342 {
343  f1 = f1 & f2;
344  return f1;
345 }
346 
347 
348 
359 namespace CellSimilarity
360 {
362  {
380  };
381 }
382 
383 
384 namespace internal
385 {
386  namespace FEValues
387  {
400  template <int dim, int spacedim=dim>
402  {
403  public:
407  void initialize (const unsigned int n_quadrature_points,
408  const UpdateFlags flags);
409 
414  std::size_t memory_consumption () const;
415 
429  std::vector<double> JxW_values;
430 
434  std::vector< DerivativeForm<1,dim,spacedim> > jacobians;
435 
440  std::vector<DerivativeForm<2,dim,spacedim> > jacobian_grads;
441 
445  std::vector<DerivativeForm<1,spacedim,dim> > inverse_jacobians;
446 
451  std::vector<Tensor<3,spacedim> > jacobian_pushed_forward_grads;
452 
457  std::vector<DerivativeForm<3,dim,spacedim> > jacobian_2nd_derivatives;
458 
463  std::vector<Tensor<4,spacedim> > jacobian_pushed_forward_2nd_derivatives;
464 
469  std::vector<DerivativeForm<4,dim,spacedim> > jacobian_3rd_derivatives;
470 
475  std::vector<Tensor<5,spacedim> > jacobian_pushed_forward_3rd_derivatives;
476 
482  std::vector<Point<spacedim> > quadrature_points;
483 
487  std::vector<Tensor<1,spacedim> > normal_vectors;
488 
492  std::vector<Tensor<1,spacedim> > boundary_forms;
493  };
494 
495 
504  template <int dim, int spacedim=dim>
506  {
507  public:
511  void initialize (const unsigned int n_quadrature_points,
512  const FiniteElement<dim,spacedim> &fe,
513  const UpdateFlags flags);
514 
519  std::size_t memory_consumption () const;
520 
539  typedef ::Table<2,double> ShapeVector;
540 
545  typedef ::Table<2,Tensor<1,spacedim> > GradientVector;
546 
550  typedef ::Table<2,Tensor<2,spacedim> > HessianVector;
551 
555  typedef ::Table<2,Tensor<3,spacedim> > ThirdDerivativeVector;
556 
562  ShapeVector shape_values;
563 
569  GradientVector shape_gradients;
570 
576  HessianVector shape_hessians;
577 
583  ThirdDerivativeVector shape_3rd_derivatives;
584 
614  std::vector<unsigned int> shape_function_to_row_table;
615  };
616  }
617 }
618 
619 
624 DEAL_II_NAMESPACE_CLOSE
625 
626 #endif
Transformed quadrature weights.
Shape function values.
Contravariant transformation.
Volume element.
Outer normal vector, not normalized.
Determinant of the Jacobian.
Transformed quadrature points.
std::vector< unsigned int > shape_function_to_row_table
void initialize(const unsigned int n_quadrature_points, const FiniteElement< dim, spacedim > &fe, const UpdateFlags flags)
Definition: fe_values.cc:2606
std::vector< Tensor< 1, spacedim > > boundary_forms
std::vector< Tensor< 4, spacedim > > jacobian_pushed_forward_2nd_derivatives
::Table< 2, Tensor< 2, spacedim > > HessianVector
std::vector< DerivativeForm< 2, dim, spacedim > > jacobian_grads
UpdateFlags & operator&=(UpdateFlags &f1, const UpdateFlags f2)
UpdateFlags operator&(const UpdateFlags f1, const UpdateFlags f2)
UpdateFlags operator|(const UpdateFlags f1, const UpdateFlags f2)
Shape function gradients of transformation.
std::vector< Tensor< 3, spacedim > > jacobian_pushed_forward_grads
No update.
::Table< 2, Tensor< 3, spacedim > > ThirdDerivativeVector
Third derivatives of shape functions.
UpdateFlags
::Table< 2, Tensor< 1, spacedim > > GradientVector
void initialize(const unsigned int n_quadrature_points, const UpdateFlags flags)
Definition: fe_values.cc:2530
std::vector< DerivativeForm< 1, dim, spacedim > > jacobians
std::vector< Tensor< 1, spacedim > > normal_vectors
Shape function values of transformation.
Second derivatives of shape functions.
Gradient of volume element.
std::vector< Tensor< 5, spacedim > > jacobian_pushed_forward_3rd_derivatives
OutputOperator< VectorType > & operator<<(OutputOperator< VectorType > &out, unsigned int step)
Definition: operator.h:163
UpdateFlags & operator|=(UpdateFlags &f1, const UpdateFlags f2)
std::vector< DerivativeForm< 4, dim, spacedim > > jacobian_3rd_derivatives
Shape function gradients.
std::vector< Point< spacedim > > quadrature_points
Normal vectors.
std::size_t memory_consumption() const
Definition: fe_values.cc:2585
Definition: fe.h:33
std::vector< DerivativeForm< 1, spacedim, dim > > inverse_jacobians
std::vector< DerivativeForm< 3, dim, spacedim > > jacobian_2nd_derivatives
Values needed for Piola transform.
Covariant transformation.