Reference documentation for deal.II version 9.5.0
\(\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\}}\)
Loading...
Searching...
No Matches
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
34template <int dim, int spacedim>
35class FiniteElement;
36
37
38namespace 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
153
161
169
177
183
189
196
203
210
218
226
242
247
253 std::array<AlignedVector<Number>, 2> shape_data_on_face;
254
264 std::array<AlignedVector<Number>, 2> quadrature_data_on_face;
265
270 std::array<AlignedVector<Number>, 2> values_within_subface;
271
276 std::array<AlignedVector<Number>, 2> gradients_within_subface;
277
282 std::array<AlignedVector<Number>, 2> hessians_within_subface;
283
288 std::array<AlignedVector<Number>, 2> subface_interpolation_matrices;
289
294 std::array<AlignedVector<typename ::internal::VectorizedArrayTrait<
295 Number>::value_type>,
296 2>
298
304
308 unsigned int fe_degree;
309
313 unsigned int n_q_points_1d;
314
320
327
334 };
335
336
337
348 template <typename Number>
350 {
357
362
366 template <int dim, int spacedim, int dim_q>
369 const unsigned int base_element = 0);
370
379 template <int dim, int spacedim, int dim_q>
380 void
382 const FiniteElement<dim, spacedim> &fe_dim,
383 const unsigned int base_element = 0);
384
401 template <int dim, int spacedim>
402 static bool
404
410 compute_orientation_table(const unsigned int n_points_per_dim);
411
419 get_shape_data(const unsigned int dimension = 0,
420 const unsigned int component = 0) const;
421
425 std::size_t
427
435 std::vector<unsigned int> lexicographic_numbering;
436
441 std::vector<UnivariateShapeData<Number>> data;
442
449
453 unsigned int n_dimensions;
454
459 unsigned int n_components;
460
465 unsigned int n_q_points;
466
472
476 unsigned int n_q_points_face;
477
482 std::vector<unsigned int> n_q_points_faces;
483
488
525
565
573
581
582 private:
587 bool
589 UnivariateShapeData<Number> &univariate_shape_data);
590
597 bool
599 const UnivariateShapeData<Number> &univariate_shape_data) const;
600 };
601
602
603
604 // ------------------------------------------ inline functions
605
606 template <typename Number>
607 inline const UnivariateShapeData<Number> &
608 ShapeInfo<Number>::get_shape_data(const unsigned int dimension,
609 const unsigned int component) const
610 {
611 AssertDimension(n_dimensions, data_access.size(0));
612 AssertDimension(n_components, data_access.size(1));
613 AssertIndexRange(dimension, n_dimensions);
614 AssertIndexRange(component, n_components);
615 return *(data_access(dimension, component));
616 }
617
618 } // end of namespace MatrixFreeFunctions
619
620} // end of namespace internal
621
623
624#endif
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static Table< 2, unsigned int > compute_orientation_table(const unsigned int n_points_per_dim)
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:524
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:572
bool check_1d_shapes_collocation(const UnivariateShapeData< Number > &univariate_shape_data) const
std::vector< UnivariateShapeData< Number > > data
Definition shape_info.h:441
::Table< 2, unsigned int > face_orientations_quad
Definition shape_info.h:580
const UnivariateShapeData< Number > & get_shape_data(const unsigned int dimension=0, const unsigned int component=0) const
Definition shape_info.h:608
std::vector< unsigned int > lexicographic_numbering
Definition shape_info.h:435
::Table< 2, unsigned int > face_to_cell_index_hermite
Definition shape_info.h:564
std::vector< unsigned int > n_q_points_faces
Definition shape_info.h:482
::Table< 2, UnivariateShapeData< Number > * > data_access
Definition shape_info.h:448
bool check_1d_shapes_symmetric(UnivariateShapeData< Number > &univariate_shape_data)
std::array< AlignedVector< Number >, 2 > shape_data_on_face
Definition shape_info.h:253
std::array< AlignedVector< typename ::internal::VectorizedArrayTrait< Number >::value_type >, 2 > subface_interpolation_matrices_scalar
Definition shape_info.h:297
std::array< AlignedVector< Number >, 2 > subface_interpolation_matrices
Definition shape_info.h:288
std::array< AlignedVector< Number >, 2 > hessians_within_subface
Definition shape_info.h:282
std::array< AlignedVector< Number >, 2 > values_within_subface
Definition shape_info.h:270
std::array< AlignedVector< Number >, 2 > quadrature_data_on_face
Definition shape_info.h:264
std::array< AlignedVector< Number >, 2 > gradients_within_subface
Definition shape_info.h:276