Reference documentation for deal.II version GIT 20f059c89a 2022-12-04 14:55: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\}}\)
qprojector.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2005 - 2021 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
30 class ReferenceCell;
31 #endif
32 
81 template <int dim>
83 {
84 public:
89  using SubQuadrature = Quadrature<dim - 1>;
90 
96  static void
98  const SubQuadrature & quadrature,
99  const unsigned int face_no,
100  std::vector<Point<dim>> &q_points);
101 
107  static Quadrature<dim>
109  const SubQuadrature &quadrature,
110  const unsigned int face_no);
111 
121  static void
123  const SubQuadrature & quadrature,
124  const unsigned int face_no,
125  const unsigned int subface_no,
126  std::vector<Point<dim>> & q_points,
127  const RefinementCase<dim - 1> &ref_case =
129 
143  static Quadrature<dim>
145  const SubQuadrature & quadrature,
146  const unsigned int face_no,
147  const unsigned int subface_no,
148  const RefinementCase<dim - 1> &ref_case =
150 
167  static Quadrature<dim>
169  const hp::QCollection<dim - 1> &quadrature);
170 
175  static Quadrature<dim>
177  const Quadrature<dim - 1> &quadrature);
178 
192  static Quadrature<dim>
194  const SubQuadrature &quadrature);
195 
206  static Quadrature<dim>
208  const Quadrature<dim> &quadrature,
209  const unsigned int child_no);
210 
220  static Quadrature<dim>
222  const Quadrature<dim> &quadrature);
223 
228  static Quadrature<dim>
230  const Quadrature<1> &quadrature,
231  const Point<dim> & p1,
232  const Point<dim> & p2);
233 
245  {
246  public:
253 
260  static DataSetDescriptor
261  cell();
262 
273  static DataSetDescriptor
275  const unsigned int face_no,
276  const bool face_orientation,
277  const bool face_flip,
278  const bool face_rotation,
279  const unsigned int n_quadrature_points);
280 
285  static DataSetDescriptor
287  const unsigned int face_no,
288  const bool face_orientation,
289  const bool face_flip,
290  const bool face_rotation,
291  const hp::QCollection<dim - 1> &quadrature);
292 
305  static DataSetDescriptor
307  const unsigned int face_no,
308  const unsigned int subface_no,
309  const bool face_orientation,
310  const bool face_flip,
311  const bool face_rotation,
312  const unsigned int n_quadrature_points,
313  const internal::SubfaceCase<dim> ref_case =
315 
322  operator unsigned int() const;
323 
324  private:
328  const unsigned int dataset_offset;
329 
334  DataSetDescriptor(const unsigned int dataset_offset);
335  };
336 };
337 
341 // ------------------- inline and template functions ----------------
342 
343 
344 
345 template <int dim>
347  const unsigned int dataset_offset)
348  : dataset_offset(dataset_offset)
349 {}
350 
351 
352 template <int dim>
354  : dataset_offset(numbers::invalid_unsigned_int)
355 {}
356 
357 
358 
359 template <int dim>
362 {
363  return 0;
364 }
365 
366 
367 
368 template <int dim>
370 {
371  return dataset_offset;
372 }
373 
374 
375 
376 template <int dim>
379  const Quadrature<dim - 1> &quadrature)
380 {
382  hp::QCollection<dim - 1>(quadrature));
383 }
384 
385 
386 /* -------------- declaration of explicit specializations ------------- */
387 
388 #ifndef DOXYGEN
389 
390 
391 template <>
392 void
394  const Quadrature<0> &,
395  const unsigned int,
396  std::vector<Point<1>> &);
397 template <>
398 void
400  const Quadrature<1> & quadrature,
401  const unsigned int face_no,
402  std::vector<Point<2>> &q_points);
403 template <>
404 void
406  const Quadrature<2> & quadrature,
407  const unsigned int face_no,
408  std::vector<Point<3>> &q_points);
409 
410 template <>
413  const hp::QCollection<0> &quadrature);
414 
415 
416 template <>
417 void
419  const Quadrature<0> &,
420  const unsigned int,
421  const unsigned int,
422  std::vector<Point<1>> &,
423  const RefinementCase<0> &);
424 template <>
425 void
427  const Quadrature<1> & quadrature,
428  const unsigned int face_no,
429  const unsigned int subface_no,
430  std::vector<Point<2>> &q_points,
431  const RefinementCase<1> &);
432 template <>
433 void
435  const Quadrature<2> & quadrature,
436  const unsigned int face_no,
437  const unsigned int subface_no,
438  std::vector<Point<3>> & q_points,
439  const RefinementCase<2> &face_ref_case);
440 
441 template <>
444  const Quadrature<0> &quadrature);
445 
446 
447 #endif // DOXYGEN
449 
450 #endif
Definition: point.h:111
const unsigned int dataset_offset
Definition: qprojector.h:328
static DataSetDescriptor cell()
Definition: qprojector.h:361
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)
Definition: qprojector.cc:1224
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_child(const ReferenceCell &reference_cell, const Quadrature< dim > &quadrature, const unsigned int child_no)
Definition: qprojector.cc:1132
static Quadrature< dim > project_to_line(const ReferenceCell &reference_cell, const Quadrature< 1 > &quadrature, const Point< dim > &p1, const Point< dim > &p2)
Definition: qprojector.cc:1196
static Quadrature< dim > project_to_all_subfaces(const ReferenceCell &reference_cell, const SubQuadrature &quadrature)
static Quadrature< dim > project_to_all_faces(const ReferenceCell &reference_cell, const hp::QCollection< dim - 1 > &quadrature)
static Quadrature< dim > project_to_all_children(const ReferenceCell &reference_cell, const Quadrature< dim > &quadrature)
Definition: qprojector.cc:1164
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:458
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:459
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
static const unsigned int invalid_unsigned_int
Definition: types.h:206