Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
shape_info.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2011 - 2019 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/aligned_vector.h>
22 #include <deal.II/base/exceptions.h>
23 #include <deal.II/base/quadrature_lib.h>
24 
25 #include <deal.II/fe/fe.h>
26 
27 
28 DEAL_II_NAMESPACE_OPEN
29 
30 namespace internal
31 {
32  namespace MatrixFreeFunctions
33  {
51  {
61 
68 
75 
80 
86 
93  };
94 
105  template <typename Number>
106  struct ShapeInfo
107  {
111  ShapeInfo();
112 
116  template <int dim>
117  ShapeInfo(const Quadrature<1> & quad,
118  const FiniteElement<dim> &fe,
119  const unsigned int base_element = 0);
120 
129  template <int dim>
130  void
131  reinit(const Quadrature<1> & quad,
132  const FiniteElement<dim> &fe_dim,
133  const unsigned int base_element = 0);
134 
138  std::size_t
139  memory_consumption() const;
140 
147 
156 
165 
174 
181 
188 
195 
202 
209 
216 
222 
228 
234 
242  std::vector<unsigned int> lexicographic_numbering;
243 
247  unsigned int fe_degree;
248 
252  unsigned int n_q_points_1d;
253 
258  unsigned int n_q_points;
259 
265 
269  unsigned int n_q_points_face;
270 
275 
281 
318 
358 
363  bool
364  check_1d_shapes_symmetric(const unsigned int n_q_points_1d);
365 
372  bool
374  };
375 
376 
377 
378  // ------------------------------------------ inline functions
379 
380  template <typename Number>
381  template <int dim>
383  const FiniteElement<dim> &fe_in,
384  const unsigned int base_element_number)
385  : element_type(tensor_general)
386  , fe_degree(0)
387  , n_q_points_1d(0)
388  , n_q_points(0)
389  , dofs_per_component_on_cell(0)
390  , n_q_points_face(0)
391  , dofs_per_component_on_face(0)
392  , nodal_at_cell_boundaries(false)
393  {
394  reinit(quad, fe_in, base_element_number);
395  }
396 
397  } // end of namespace MatrixFreeFunctions
398 
399 } // end of namespace internal
400 
401 DEAL_II_NAMESPACE_CLOSE
402 
403 #endif
void reinit(const Quadrature< 1 > &quad, const FiniteElement< dim > &fe_dim, const unsigned int base_element=0)
AlignedVector< Number > shape_values_eo
Definition: shape_info.h:180
AlignedVector< Number > shape_hessians
Definition: shape_info.h:173
AlignedVector< Number > shape_values
Definition: shape_info.h:155
AlignedVector< Number > hessians_within_subface[2]
Definition: shape_info.h:233
bool check_1d_shapes_symmetric(const unsigned int n_q_points_1d)
AlignedVector< Number > shape_gradients_collocation_eo
Definition: shape_info.h:201
AlignedVector< Number > shape_data_on_face[2]
Definition: shape_info.h:215
AlignedVector< Number > shape_hessians_eo
Definition: shape_info.h:194
::Table< 2, unsigned int > face_to_cell_index_hermite
Definition: shape_info.h:357
AlignedVector< Number > gradients_within_subface[2]
Definition: shape_info.h:227
::Table< 2, unsigned int > face_to_cell_index_nodal
Definition: shape_info.h:317
AlignedVector< Number > shape_gradients_eo
Definition: shape_info.h:187
AlignedVector< Number > shape_gradients
Definition: shape_info.h:164
AlignedVector< Number > values_within_subface[2]
Definition: shape_info.h:221
std::vector< unsigned int > lexicographic_numbering
Definition: shape_info.h:242
AlignedVector< Number > shape_hessians_collocation_eo
Definition: shape_info.h:208