Reference documentation for deal.II version GIT relicensing-136-gb80d0be4af 2024-03-18 08:20:02+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\}}\)
Loading...
Searching...
No Matches
qprojector.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2005 - 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_qprojector_h
16#define dealii_qprojector_h
17
18
19#include <deal.II/base/config.h>
20
23
25
27
28#ifndef DOXYGEN
29class ReferenceCell;
30#endif
31
80template <int dim>
82{
83public:
88 using SubQuadrature = Quadrature<dim - 1>;
89
95 static void
96 project_to_face(const ReferenceCell &reference_cell,
97 const SubQuadrature &quadrature,
98 const unsigned int face_no,
99 std::vector<Point<dim>> &q_points);
100
106 static Quadrature<dim>
107 project_to_face(const ReferenceCell &reference_cell,
108 const SubQuadrature &quadrature,
109 const unsigned int face_no);
110
117 static Quadrature<dim>
118 project_to_oriented_face(const ReferenceCell &reference_cell,
119 const SubQuadrature &quadrature,
120 const unsigned int face_no,
121 const bool face_orientation,
122 const bool face_flip,
123 const bool face_rotation);
124
134 static void
135 project_to_subface(const ReferenceCell &reference_cell,
136 const SubQuadrature &quadrature,
137 const unsigned int face_no,
138 const unsigned int subface_no,
139 std::vector<Point<dim>> &q_points,
140 const RefinementCase<dim - 1> &ref_case =
142
156 static Quadrature<dim>
157 project_to_subface(const ReferenceCell &reference_cell,
158 const SubQuadrature &quadrature,
159 const unsigned int face_no,
160 const unsigned int subface_no,
161 const RefinementCase<dim - 1> &ref_case =
163
173 static Quadrature<dim>
174 project_to_oriented_subface(const ReferenceCell &reference_cell,
175 const SubQuadrature &quadrature,
176 const unsigned int face_no,
177 const unsigned int subface_no,
178 const bool face_orientation,
179 const bool face_flip,
180 const bool face_rotation,
181 const internal::SubfaceCase<dim> ref_case);
182
199 static Quadrature<dim>
200 project_to_all_faces(const ReferenceCell &reference_cell,
201 const hp::QCollection<dim - 1> &quadrature);
202
207 static Quadrature<dim>
208 project_to_all_faces(const ReferenceCell &reference_cell,
209 const Quadrature<dim - 1> &quadrature);
210
224 static Quadrature<dim>
226 const SubQuadrature &quadrature);
227
238 static Quadrature<dim>
239 project_to_child(const ReferenceCell &reference_cell,
240 const Quadrature<dim> &quadrature,
241 const unsigned int child_no);
242
252 static Quadrature<dim>
253 project_to_all_children(const ReferenceCell &reference_cell,
254 const Quadrature<dim> &quadrature);
255
260 static Quadrature<dim>
261 project_to_line(const ReferenceCell &reference_cell,
262 const Quadrature<1> &quadrature,
263 const Point<dim> &p1,
264 const Point<dim> &p2);
265
277 {
278 public:
285
292 static DataSetDescriptor
293 cell();
294
305 static DataSetDescriptor
306 face(const ReferenceCell &reference_cell,
307 const unsigned int face_no,
308 const bool face_orientation,
309 const bool face_flip,
310 const bool face_rotation,
311 const unsigned int n_quadrature_points);
312
317 static DataSetDescriptor
318 face(const ReferenceCell &reference_cell,
319 const unsigned int face_no,
320 const bool face_orientation,
321 const bool face_flip,
322 const bool face_rotation,
323 const hp::QCollection<dim - 1> &quadrature);
324
337 static DataSetDescriptor
338 subface(const ReferenceCell &reference_cell,
339 const unsigned int face_no,
340 const unsigned int subface_no,
341 const bool face_orientation,
342 const bool face_flip,
343 const bool face_rotation,
344 const unsigned int n_quadrature_points,
345 const internal::SubfaceCase<dim> ref_case =
347
354 operator unsigned int() const;
355
356 private:
360 const unsigned int dataset_offset;
361
366 DataSetDescriptor(const unsigned int dataset_offset);
367 };
368};
369
373// ------------------- inline and template functions ----------------
374
375
376
377template <int dim>
379 const unsigned int dataset_offset)
380 : dataset_offset(dataset_offset)
381{}
382
383
384template <int dim>
386 : dataset_offset(numbers::invalid_unsigned_int)
387{}
388
389
390
391template <int dim>
397
398
399
400template <int dim>
402{
403 return dataset_offset;
404}
405
406
407
408template <int dim>
410 const ReferenceCell &reference_cell,
411 const Quadrature<dim - 1> &quadrature)
412{
413 return project_to_all_faces(reference_cell,
414 hp::QCollection<dim - 1>(quadrature));
415}
416
417
418/* -------------- declaration of explicit specializations ------------- */
419
420#ifndef DOXYGEN
421
422
423template <>
424void
426 const Quadrature<0> &,
427 const unsigned int,
428 std::vector<Point<1>> &);
429template <>
430void
432 const Quadrature<1> &quadrature,
433 const unsigned int face_no,
434 std::vector<Point<2>> &q_points);
435template <>
436void
438 const Quadrature<2> &quadrature,
439 const unsigned int face_no,
440 std::vector<Point<3>> &q_points);
441
442template <>
445 const hp::QCollection<0> &quadrature);
446
447
448template <>
449void
451 const Quadrature<0> &,
452 const unsigned int,
453 const unsigned int,
454 std::vector<Point<1>> &,
455 const RefinementCase<0> &);
456template <>
457void
459 const Quadrature<1> &quadrature,
460 const unsigned int face_no,
461 const unsigned int subface_no,
462 std::vector<Point<2>> &q_points,
463 const RefinementCase<1> &);
464template <>
465void
467 const Quadrature<2> &quadrature,
468 const unsigned int face_no,
469 const unsigned int subface_no,
470 std::vector<Point<3>> &q_points,
471 const RefinementCase<2> &face_ref_case);
472
473template <>
476 const Quadrature<0> &quadrature);
477
478
479#endif // DOXYGEN
481
482#endif
Definition point.h:111
const unsigned int dataset_offset
Definition qprojector.h:360
static DataSetDescriptor cell()
Definition qprojector.h:393
static DataSetDescriptor subface(const ReferenceCell &reference_cell, const unsigned int face_no, const unsigned int subface_no, const bool face_orientation, const bool face_flip, const bool face_rotation, const unsigned int n_quadrature_points, const internal::SubfaceCase< dim > ref_case=internal::SubfaceCase< dim >::case_isotropic)
static DataSetDescriptor face(const ReferenceCell &reference_cell, const unsigned int face_no, const bool face_orientation, const bool face_flip, const bool face_rotation, const unsigned int n_quadrature_points)
static Quadrature< dim > project_to_all_subfaces(const ReferenceCell &reference_cell, const SubQuadrature &quadrature)
static Quadrature< dim > project_to_child(const ReferenceCell &reference_cell, const Quadrature< dim > &quadrature, const unsigned int child_no)
static Quadrature< dim > project_to_oriented_subface(const ReferenceCell &reference_cell, const SubQuadrature &quadrature, const unsigned int face_no, const unsigned int subface_no, const bool face_orientation, const bool face_flip, const bool face_rotation, const internal::SubfaceCase< dim > ref_case)
static Quadrature< dim > project_to_all_faces(const ReferenceCell &reference_cell, const hp::QCollection< dim - 1 > &quadrature)
static Quadrature< dim > project_to_oriented_face(const ReferenceCell &reference_cell, const SubQuadrature &quadrature, const unsigned int face_no, const bool face_orientation, const bool face_flip, const bool face_rotation)
static Quadrature< dim > project_to_line(const ReferenceCell &reference_cell, const Quadrature< 1 > &quadrature, const Point< dim > &p1, const Point< dim > &p2)
static void project_to_subface(const ReferenceCell &reference_cell, const SubQuadrature &quadrature, const unsigned int face_no, const unsigned int subface_no, std::vector< Point< dim > > &q_points, const RefinementCase< dim - 1 > &ref_case=RefinementCase< dim - 1 >::isotropic_refinement)
static Quadrature< dim > project_to_all_children(const ReferenceCell &reference_cell, const Quadrature< dim > &quadrature)
static void project_to_face(const ReferenceCell &reference_cell, const SubQuadrature &quadrature, const unsigned int face_no, std::vector< Point< dim > > &q_points)
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503