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
util.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2020 - 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_util_h
18#define dealii_matrix_free_util_h
19
20
21#include <deal.II/base/config.h>
22
25
27
29
31
32
33namespace internal
34{
35 namespace MatrixFreeFunctions
36 {
42 template <int dim>
43 inline std::pair<ReferenceCell, ::hp::QCollection<dim - 1>>
45 const bool do_assert = true)
46 {
47 if (dim == 2 || dim == 3)
48 {
49 for (unsigned int i = 1; i <= 4; ++i)
50 if (quad == QGaussSimplex<dim>(i))
51 return {ReferenceCells::get_simplex<dim>(),
52 ::hp::QCollection<dim - 1>(
53 QGaussSimplex<dim - 1>(i))};
54
55 for (unsigned int i = 1; i <= 5; ++i)
56 if (quad == QWitherdenVincentSimplex<dim>(i))
57 return {ReferenceCells::get_simplex<dim>(),
58 ::hp::QCollection<dim - 1>(
59 QWitherdenVincentSimplex<dim - 1>(i))};
60 }
61
62 if (dim == 3)
63 for (unsigned int i = 1; i <= 3; ++i)
64 if (quad == QGaussWedge<dim>(i))
65 {
66 QGauss<dim - 1> quad(i);
67 QGaussSimplex<dim - 1> tri(i);
68
69 return {
71 ::hp::QCollection<dim - 1>(tri, tri, quad, quad, quad)};
72 }
73
74 if (dim == 3)
75 for (unsigned int i = 1; i <= 2; ++i)
76 if (quad == QGaussPyramid<dim>(i))
77 {
78 QGauss<dim - 1> quad(i);
79 QGaussSimplex<dim - 1> tri(i);
80
81 return {
83 ::hp::QCollection<dim - 1>(quad, tri, tri, tri, tri)};
84 }
85
86 // note: handle hypercubes last since normally this function is not
87 // called for hypercubes
88 for (unsigned int i = 1; i <= 5; ++i)
89 if (quad == QGauss<dim>(i))
90 return {ReferenceCells::get_hypercube<dim>(),
91 ::hp::QCollection<dim - 1>(QGauss<dim - 1>(i))};
92
93 if (do_assert)
95
96 return {ReferenceCells::Invalid, ::hp::QCollection<dim - 1>()};
97 }
98
99
100
109 template <int dim>
110 inline std::pair<Quadrature<dim - 1>, Quadrature<dim - 1>>
112 {
114 quad.size() > 0,
116 "There is nothing useful you can do with a MatrixFree/FEEvaluation "
117 "object when using a quadrature formula with zero "
118 "quadrature points!"));
119
120 if (dim == 2 || dim == 3)
121 {
122 for (unsigned int i = 1; i <= 4; ++i)
123 if (quad == QGaussSimplex<dim>(i))
124 {
125 if (dim == 2)
126 return {QGaussSimplex<dim - 1>(i), // line!
128 else
129 return {Quadrature<dim - 1>(), QGaussSimplex<dim - 1>(i)};
130 }
131
132 for (unsigned int i = 1; i <= 5; ++i)
133 if (quad == QWitherdenVincentSimplex<dim>(i))
134 {
135 if (dim == 2)
136 return {QWitherdenVincentSimplex<dim - 1>(i), // line!
138 else
139 return {Quadrature<dim - 1>(),
141 }
142 }
143
144 if (dim == 3)
145 for (unsigned int i = 1; i <= 3; ++i)
146 if (quad == QGaussWedge<dim>(i))
147 return {QGauss<dim - 1>(i), QGaussSimplex<dim - 1>(i)};
148
149 if (dim == 3)
150 for (unsigned int i = 1; i <= 2; ++i)
151 if (quad == QGaussPyramid<dim>(i))
152 return {QGauss<dim - 1>(i), QGaussSimplex<dim - 1>(i)};
153
154 // note: handle hypercubes last since normally this function is not
155 // called for hypercubes
156 for (unsigned int i = 1; i <= 5; ++i)
157 if (quad == QGauss<dim>(i))
158 return {QGauss<dim - 1>(i), Quadrature<dim - 1>()};
159
161
162 return {Quadrature<dim - 1>(), Quadrature<dim - 1>()};
163 }
164
165 inline DEAL_II_ALWAYS_INLINE unsigned int
166 indicate_power_of_two(const unsigned int vectorization_length)
167 {
168 unsigned int vectorization_length_bits = 0;
169 unsigned int my_length = vectorization_length;
170 while (my_length >>= 1)
171 ++vectorization_length_bits;
172 return 1 << vectorization_length_bits;
173 }
174
175 } // end of namespace MatrixFreeFunctions
176} // end of namespace internal
177
179
180#endif
unsigned int size() const
#define DEAL_II_ALWAYS_INLINE
Definition config.h:106
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
static ::ExceptionBase & ExcNotImplemented()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
constexpr const ReferenceCell Wedge
constexpr const ReferenceCell Pyramid
constexpr const ReferenceCell Invalid
unsigned int indicate_power_of_two(const unsigned int vectorization_length)
Definition util.h:166
std::pair< Quadrature< dim - 1 >, Quadrature< dim - 1 > > get_unique_face_quadratures(const Quadrature< dim > &quad)
Definition util.h:111
std::pair< ReferenceCell, ::hp::QCollection< dim - 1 > > get_face_quadrature_collection(const Quadrature< dim > &quad, const bool do_assert=true)
Definition util.h:44