Reference documentation for deal.II version Git e67cfc96e3 2021-10-16 16:07:33 +0200
\(\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\}}\)
shape_info.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2011 - 2021 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 
17 #ifndef dealii_matrix_free_shape_info_h
18 #define dealii_matrix_free_shape_info_h
19 
20 
21 #include <deal.II/base/config.h>
22 
26 
27 #include <deal.II/fe/fe.h>
28 
29 
31 
32 namespace internal
33 {
34  namespace MatrixFreeFunctions
35  {
53  {
63 
70 
77 
82 
88 
95 
100  };
101 
102 
103 
112  template <typename Number>
114  {
119 
123  std::size_t
124  memory_consumption() const;
125 
132 
140 
148 
156 
162 
168 
175 
182 
189 
197 
205 
221 
226 
232  std::array<AlignedVector<Number>, 2> shape_data_on_face;
233 
244  std::array<AlignedVector<Number>, 2> quadrature_data_on_face;
245 
250  std::array<AlignedVector<Number>, 2> values_within_subface;
251 
256  std::array<AlignedVector<Number>, 2> gradients_within_subface;
257 
262  std::array<AlignedVector<Number>, 2> hessians_within_subface;
263 
268  std::array<AlignedVector<Number>, 2> subface_interpolation_matrices;
269 
275 
279  unsigned int fe_degree;
280 
284  unsigned int n_q_points_1d;
285 
291 
298 
305  };
306 
307 
308 
319  template <typename Number>
320  struct ShapeInfo
321  {
328 
332  ShapeInfo();
333 
337  template <int dim, int dim_q>
338  ShapeInfo(const Quadrature<dim_q> & quad,
339  const FiniteElement<dim> &fe,
340  const unsigned int base_element = 0);
341 
350  template <int dim, int dim_q>
351  void
352  reinit(const Quadrature<dim_q> & quad,
353  const FiniteElement<dim> &fe_dim,
354  const unsigned int base_element = 0);
355 
359  template <int dim, int spacedim>
360  static bool
361  is_supported(const FiniteElement<dim, spacedim> &fe);
362 
370  get_shape_data(const unsigned int dimension = 0,
371  const unsigned int component = 0) const;
372 
376  std::size_t
377  memory_consumption() const;
378 
386  std::vector<unsigned int> lexicographic_numbering;
387 
392  std::vector<UnivariateShapeData<Number>> data;
393 
400 
404  unsigned int n_dimensions;
405 
410  unsigned int n_components;
411 
416  unsigned int n_q_points;
417 
423 
427  unsigned int n_q_points_face;
428 
433  std::vector<unsigned int> n_q_points_faces;
434 
439 
476 
516 
524 
525  private:
530  bool
531  check_1d_shapes_symmetric(
532  UnivariateShapeData<Number> &univariate_shape_data);
533 
540  bool
541  check_1d_shapes_collocation(
542  const UnivariateShapeData<Number> &univariate_shape_data) const;
543  };
544 
545 
546 
547  // ------------------------------------------ inline functions
548 
549  template <typename Number>
550  template <int dim, int dim_q>
552  const FiniteElement<dim> &fe_in,
553  const unsigned int base_element_number)
555  , n_dimensions(0)
556  , n_components(0)
557  , n_q_points(0)
558  , dofs_per_component_on_cell(0)
559  , n_q_points_face(0)
560  , dofs_per_component_on_face(0)
561  {
562  reinit(quad, fe_in, base_element_number);
563  }
564 
565  template <typename Number>
566  inline const UnivariateShapeData<Number> &
567  ShapeInfo<Number>::get_shape_data(const unsigned int dimension,
568  const unsigned int component) const
569  {
570  AssertDimension(n_dimensions, data_access.size(0));
571  AssertDimension(n_components, data_access.size(1));
572  AssertIndexRange(dimension, n_dimensions);
573  AssertIndexRange(component, n_components);
574  return *(data_access(dimension, component));
575  }
576 
577  } // end of namespace MatrixFreeFunctions
578 
579 } // end of namespace internal
580 
582 
583 #endif
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1655
#define AssertIndexRange(index, range)
Definition: exceptions.h:1720
AlignedVector< Number > shape_hessians_collocation_eo
Definition: shape_info.h:204
std::array< AlignedVector< Number >, 2 > shape_data_on_face
Definition: shape_info.h:232
std::vector< unsigned int > n_q_points_faces
Definition: shape_info.h:433
std::array< AlignedVector< Number >, 2 > values_within_subface
Definition: shape_info.h:250
std::array< AlignedVector< Number >, 2 > gradients_within_subface
Definition: shape_info.h:256
::Table< 2, UnivariateShapeData< Number > * > data_access
Definition: shape_info.h:399
::Table< 2, unsigned int > face_orientations
Definition: shape_info.h:523
AlignedVector< Number > shape_gradients_collocation_eo
Definition: shape_info.h:196
std::array< AlignedVector< Number >, 2 > quadrature_data_on_face
Definition: shape_info.h:244
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
std::array< AlignedVector< Number >, 2 > subface_interpolation_matrices
Definition: shape_info.h:268
std::array< AlignedVector< Number >, 2 > hessians_within_subface
Definition: shape_info.h:262
::Table< 2, unsigned int > face_to_cell_index_hermite
Definition: shape_info.h:515
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
::Table< 2, unsigned int > face_to_cell_index_nodal
Definition: shape_info.h:475
std::vector< UnivariateShapeData< Number > > data
Definition: shape_info.h:392
std::vector< unsigned int > lexicographic_numbering
Definition: shape_info.h:386