Reference documentation for deal.II version GIT 29f9da0a34 2023-12-07 10:00:01+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\}}\)
shape_info.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2011 - 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 
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 #include <deal.II/base/table.h>
28 
29 
31 
32 
33 // forward declaration
34 template <int dim, int spacedim>
35 class FiniteElement;
36 
37 
38 namespace internal
39 {
40  namespace MatrixFreeFunctions
41  {
59  {
69 
76 
83 
89 
94 
100 
107 
114 
118  tensor_none = 8
119 
120 
121  };
122 
123 
124 
133  template <typename Number>
135  {
140 
144  std::size_t
146 
156  template <int dim, int spacedim>
157  void
159  const Quadrature<1> &quad,
160  const std::vector<unsigned int> &lexicographic,
161  const unsigned int direction);
162 
169  template <int dim, int spacedim>
170  void
172  const Quadrature<1> &quad,
173  const std::vector<unsigned int> &lexicographic,
174  const unsigned int direction);
175 
180  bool
182 
189  bool
191 
198 
206 
214 
222 
228 
234 
241 
248 
255 
263 
271 
287 
292 
298  std::array<AlignedVector<Number>, 2> shape_data_on_face;
299 
309  std::array<AlignedVector<Number>, 2> quadrature_data_on_face;
310 
315  std::array<AlignedVector<Number>, 2> values_within_subface;
316 
321  std::array<AlignedVector<Number>, 2> gradients_within_subface;
322 
327  std::array<AlignedVector<Number>, 2> hessians_within_subface;
328 
333  std::array<AlignedVector<Number>, 2> subface_interpolation_matrices;
334 
339  std::array<AlignedVector<typename ::internal::VectorizedArrayTrait<
340  Number>::value_type>,
341  2>
343 
349 
353  unsigned int fe_degree;
354 
358  unsigned int n_q_points_1d;
359 
365 
372 
379  };
380 
381 
382 
393  template <typename Number>
394  struct ShapeInfo
395  {
402 
407 
411  template <int dim, int spacedim, int dim_q>
414  const unsigned int base_element = 0);
415 
424  template <int dim, int spacedim, int dim_q>
425  void
427  const FiniteElement<dim, spacedim> &fe_dim,
428  const unsigned int base_element = 0);
429 
446  template <int dim, int spacedim>
447  static bool
449 
455  compute_orientation_table(const unsigned int n_points_per_dim);
456 
464  get_shape_data(const unsigned int dimension = 0,
465  const unsigned int component = 0) const;
466 
470  std::size_t
472 
480  std::vector<unsigned int> lexicographic_numbering;
481 
486  std::vector<UnivariateShapeData<Number>> data;
487 
494 
498  unsigned int n_dimensions;
499 
504  unsigned int n_components;
505 
510  unsigned int n_q_points;
511 
517 
521  unsigned int n_q_points_face;
522 
527  std::vector<unsigned int> n_q_points_faces;
528 
533 
570 
610 
618 
626  };
627 
628 
629 
630  // ------------------------------------------ inline functions
631 
632  template <typename Number>
633  inline const UnivariateShapeData<Number> &
634  ShapeInfo<Number>::get_shape_data(const unsigned int dimension,
635  const unsigned int component) const
636  {
637  AssertDimension(n_dimensions, data_access.size(0));
638  AssertDimension(n_components, data_access.size(1));
639  AssertIndexRange(dimension, n_dimensions);
640  AssertIndexRange(component, n_components);
641  return *(data_access(dimension, component));
642  }
643 
644  } // end of namespace MatrixFreeFunctions
645 
646 } // end of namespace internal
647 
649 
650 #endif
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1820
#define AssertIndexRange(index, range)
Definition: exceptions.h:1888
ShapeInfo(const Quadrature< dim_q > &quad, const FiniteElement< dim, spacedim > &fe, const unsigned int base_element=0)
::Table< 2, unsigned int > face_to_cell_index_nodal
Definition: shape_info.h:569
void reinit(const Quadrature< dim_q > &quad, const FiniteElement< dim, spacedim > &fe_dim, const unsigned int base_element=0)
static bool is_supported(const FiniteElement< dim, spacedim > &fe)
::Table< 2, unsigned int > face_orientations_dofs
Definition: shape_info.h:617
std::vector< UnivariateShapeData< Number > > data
Definition: shape_info.h:486
::Table< 2, unsigned int > face_orientations_quad
Definition: shape_info.h:625
const UnivariateShapeData< Number > & get_shape_data(const unsigned int dimension=0, const unsigned int component=0) const
Definition: shape_info.h:634
std::vector< unsigned int > lexicographic_numbering
Definition: shape_info.h:480
::Table< 2, unsigned int > face_to_cell_index_hermite
Definition: shape_info.h:609
std::vector< unsigned int > n_q_points_faces
Definition: shape_info.h:527
::Table< 2, UnivariateShapeData< Number > * > data_access
Definition: shape_info.h:493
static Table< 2, unsigned int > compute_orientation_table(const unsigned int n_points_per_dim)
void evaluate_collocation_space(const FiniteElement< dim, spacedim > &fe, const Quadrature< 1 > &quad, const std::vector< unsigned int > &lexicographic, const unsigned int direction)
std::array< AlignedVector< Number >, 2 > shape_data_on_face
Definition: shape_info.h:298
std::array< AlignedVector< typename ::internal::VectorizedArrayTrait< Number >::value_type >, 2 > subface_interpolation_matrices_scalar
Definition: shape_info.h:342
AlignedVector< Number > shape_gradients_collocation_eo
Definition: shape_info.h:262
std::array< AlignedVector< Number >, 2 > subface_interpolation_matrices
Definition: shape_info.h:333
void evaluate_shape_functions(const FiniteElement< dim, spacedim > &fe, const Quadrature< 1 > &quad, const std::vector< unsigned int > &lexicographic, const unsigned int direction)
std::array< AlignedVector< Number >, 2 > hessians_within_subface
Definition: shape_info.h:327
std::array< AlignedVector< Number >, 2 > values_within_subface
Definition: shape_info.h:315
std::array< AlignedVector< Number >, 2 > quadrature_data_on_face
Definition: shape_info.h:309
AlignedVector< Number > shape_hessians_collocation_eo
Definition: shape_info.h:270
std::array< AlignedVector< Number >, 2 > gradients_within_subface
Definition: shape_info.h:321