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
qprojector.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2005 - 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#ifndef dealii_qprojector_h
17#define dealii_qprojector_h
18
19
20#include <deal.II/base/config.h>
21
24
26
28
29#ifndef DOXYGEN
30class ReferenceCell;
31#endif
32
81template <int dim>
83{
84public:
89 using SubQuadrature = Quadrature<dim - 1>;
90
96 static void
97 project_to_face(const ReferenceCell & reference_cell,
98 const SubQuadrature & quadrature,
99 const unsigned int face_no,
100 std::vector<Point<dim>> &q_points);
101
107 static Quadrature<dim>
108 project_to_face(const ReferenceCell &reference_cell,
109 const SubQuadrature &quadrature,
110 const unsigned int face_no);
111
118 static Quadrature<dim>
119 project_to_oriented_face(const ReferenceCell &reference_cell,
120 const SubQuadrature &quadrature,
121 const unsigned int face_no,
122 const bool face_orientation,
123 const bool face_flip,
124 const bool face_rotation);
125
135 static void
136 project_to_subface(const ReferenceCell & reference_cell,
137 const SubQuadrature & quadrature,
138 const unsigned int face_no,
139 const unsigned int subface_no,
140 std::vector<Point<dim>> & q_points,
141 const RefinementCase<dim - 1> &ref_case =
143
157 static Quadrature<dim>
158 project_to_subface(const ReferenceCell & reference_cell,
159 const SubQuadrature & quadrature,
160 const unsigned int face_no,
161 const unsigned int subface_no,
162 const RefinementCase<dim - 1> &ref_case =
164
174 static Quadrature<dim>
175 project_to_oriented_subface(const ReferenceCell & reference_cell,
176 const SubQuadrature & quadrature,
177 const unsigned int face_no,
178 const unsigned int subface_no,
179 const bool face_orientation,
180 const bool face_flip,
181 const bool face_rotation,
182 const internal::SubfaceCase<dim> ref_case);
183
200 static Quadrature<dim>
201 project_to_all_faces(const ReferenceCell & reference_cell,
202 const hp::QCollection<dim - 1> &quadrature);
203
208 static Quadrature<dim>
209 project_to_all_faces(const ReferenceCell & reference_cell,
210 const Quadrature<dim - 1> &quadrature);
211
225 static Quadrature<dim>
227 const SubQuadrature &quadrature);
228
239 static Quadrature<dim>
240 project_to_child(const ReferenceCell & reference_cell,
241 const Quadrature<dim> &quadrature,
242 const unsigned int child_no);
243
253 static Quadrature<dim>
254 project_to_all_children(const ReferenceCell & reference_cell,
255 const Quadrature<dim> &quadrature);
256
261 static Quadrature<dim>
262 project_to_line(const ReferenceCell &reference_cell,
263 const Quadrature<1> &quadrature,
264 const Point<dim> & p1,
265 const Point<dim> & p2);
266
278 {
279 public:
286
293 static DataSetDescriptor
294 cell();
295
306 static DataSetDescriptor
307 face(const ReferenceCell &reference_cell,
308 const unsigned int face_no,
309 const bool face_orientation,
310 const bool face_flip,
311 const bool face_rotation,
312 const unsigned int n_quadrature_points);
313
318 static DataSetDescriptor
319 face(const ReferenceCell & reference_cell,
320 const unsigned int face_no,
321 const bool face_orientation,
322 const bool face_flip,
323 const bool face_rotation,
324 const hp::QCollection<dim - 1> &quadrature);
325
338 static DataSetDescriptor
339 subface(const ReferenceCell & reference_cell,
340 const unsigned int face_no,
341 const unsigned int subface_no,
342 const bool face_orientation,
343 const bool face_flip,
344 const bool face_rotation,
345 const unsigned int n_quadrature_points,
346 const internal::SubfaceCase<dim> ref_case =
348
355 operator unsigned int() const;
356
357 private:
361 const unsigned int dataset_offset;
362
367 DataSetDescriptor(const unsigned int dataset_offset);
368 };
369};
370
374// ------------------- inline and template functions ----------------
375
376
377
378template <int dim>
380 const unsigned int dataset_offset)
381 : dataset_offset(dataset_offset)
382{}
383
384
385template <int dim>
387 : dataset_offset(numbers::invalid_unsigned_int)
388{}
389
390
391
392template <int dim>
395{
396 return 0;
397}
398
399
400
401template <int dim>
403{
404 return dataset_offset;
405}
406
407
408
409template <int dim>
411 const ReferenceCell & reference_cell,
412 const Quadrature<dim - 1> &quadrature)
413{
414 return project_to_all_faces(reference_cell,
415 hp::QCollection<dim - 1>(quadrature));
416}
417
418
419/* -------------- declaration of explicit specializations ------------- */
420
421#ifndef DOXYGEN
422
423
424template <>
425void
427 const Quadrature<0> &,
428 const unsigned int,
429 std::vector<Point<1>> &);
430template <>
431void
432QProjector<2>::project_to_face(const ReferenceCell & reference_cell,
433 const Quadrature<1> & quadrature,
434 const unsigned int face_no,
435 std::vector<Point<2>> &q_points);
436template <>
437void
438QProjector<3>::project_to_face(const ReferenceCell & reference_cell,
439 const Quadrature<2> & quadrature,
440 const unsigned int face_no,
441 std::vector<Point<3>> &q_points);
442
443template <>
446 const hp::QCollection<0> &quadrature);
447
448
449template <>
450void
452 const Quadrature<0> &,
453 const unsigned int,
454 const unsigned int,
455 std::vector<Point<1>> &,
456 const RefinementCase<0> &);
457template <>
458void
460 const Quadrature<1> & quadrature,
461 const unsigned int face_no,
462 const unsigned int subface_no,
463 std::vector<Point<2>> &q_points,
464 const RefinementCase<1> &);
465template <>
466void
468 const Quadrature<2> & quadrature,
469 const unsigned int face_no,
470 const unsigned int subface_no,
471 std::vector<Point<3>> & q_points,
472 const RefinementCase<2> &face_ref_case);
473
474template <>
477 const Quadrature<0> &quadrature);
478
479
480#endif // DOXYGEN
482
483#endif
Definition point.h:112
const unsigned int dataset_offset
Definition qprojector.h:361
static DataSetDescriptor cell()
Definition qprojector.h:394
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:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473