Reference documentation for deal.II version Git ac854ef673 2021-03-06 10:09:45 +0100
\(\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 - 2019 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 
30 
31 
34 
35 
79 template <int dim>
81 {
82 public:
87  using SubQuadrature = Quadrature<dim - 1>;
88 
98  DEAL_II_DEPRECATED static void
99  project_to_face(const SubQuadrature & quadrature,
100  const unsigned int face_no,
101  std::vector<Point<dim>> &q_points);
102 
108  static void
110  const SubQuadrature & quadrature,
111  const unsigned int face_no,
112  std::vector<Point<dim>> &q_points);
113 
124  project_to_face(const SubQuadrature &quadrature, const unsigned int face_no);
125 
131  static Quadrature<dim>
133  const SubQuadrature &quadrature,
134  const unsigned int face_no);
135 
149  DEAL_II_DEPRECATED static void
150  project_to_subface(const SubQuadrature & quadrature,
151  const unsigned int face_no,
152  const unsigned int subface_no,
153  std::vector<Point<dim>> & q_points,
154  const RefinementCase<dim - 1> &ref_case =
156 
166  static void
168  const SubQuadrature & quadrature,
169  const unsigned int face_no,
170  const unsigned int subface_no,
171  std::vector<Point<dim>> & q_points,
172  const RefinementCase<dim - 1> &ref_case =
174 
189  project_to_subface(const SubQuadrature & quadrature,
190  const unsigned int face_no,
191  const unsigned int subface_no,
192  const RefinementCase<dim - 1> &ref_case =
194 
208  static Quadrature<dim>
210  const SubQuadrature & quadrature,
211  const unsigned int face_no,
212  const unsigned int subface_no,
213  const RefinementCase<dim - 1> &ref_case =
215 
237  project_to_all_faces(const Quadrature<dim - 1> &quadrature);
238 
255  static Quadrature<dim>
257  const hp::QCollection<dim - 1> &quadrature);
258 
263  static Quadrature<dim>
265  const Quadrature<dim - 1> &quadrature);
266 
285  project_to_all_subfaces(const SubQuadrature &quadrature);
286 
300  static Quadrature<dim>
302  const SubQuadrature &quadrature);
303 
319  project_to_child(const Quadrature<dim> &quadrature,
320  const unsigned int child_no);
321 
332  static Quadrature<dim>
334  const Quadrature<dim> &quadrature,
335  const unsigned int child_no);
336 
351  project_to_all_children(const Quadrature<dim> &quadrature);
352 
362  static Quadrature<dim>
364  const Quadrature<dim> &quadrature);
365 
375  project_to_line(const Quadrature<1> &quadrature,
376  const Point<dim> & p1,
377  const Point<dim> & p2);
378 
383  static Quadrature<dim>
385  const Quadrature<1> &quadrature,
386  const Point<dim> & p1,
387  const Point<dim> & p2);
388 
400  {
401  public:
408 
415  static DataSetDescriptor
416  cell();
417 
433  face(const unsigned int face_no,
434  const bool face_orientation,
435  const bool face_flip,
436  const bool face_rotation,
437  const unsigned int n_quadrature_points);
438 
449  static DataSetDescriptor
451  const unsigned int face_no,
452  const bool face_orientation,
453  const bool face_flip,
454  const bool face_rotation,
455  const unsigned int n_quadrature_points);
456 
461  static DataSetDescriptor
462  face(const ReferenceCell reference_cell,
463  const unsigned int face_no,
464  const bool face_orientation,
465  const bool face_flip,
466  const bool face_rotation,
467  const hp::QCollection<dim - 1> &quadrature);
468 
486  subface(const unsigned int face_no,
487  const unsigned int subface_no,
488  const bool face_orientation,
489  const bool face_flip,
490  const bool face_rotation,
491  const unsigned int n_quadrature_points,
492  const internal::SubfaceCase<dim> ref_case =
494 
507  static DataSetDescriptor
508  subface(const ReferenceCell reference_cell,
509  const unsigned int face_no,
510  const unsigned int subface_no,
511  const bool face_orientation,
512  const bool face_flip,
513  const bool face_rotation,
514  const unsigned int n_quadrature_points,
515  const internal::SubfaceCase<dim> ref_case =
517 
524  operator unsigned int() const;
525 
526  private:
530  const unsigned int dataset_offset;
531 
536  DataSetDescriptor(const unsigned int dataset_offset);
537  };
538 };
539 
543 // ------------------- inline and template functions ----------------
544 
545 
546 
547 template <int dim>
549  const unsigned int dataset_offset)
550  : dataset_offset(dataset_offset)
551 {}
552 
553 
554 template <int dim>
557 {}
558 
559 
560 
561 template <int dim>
564 {
565  return 0;
566 }
567 
568 
569 
570 template <int dim>
572 {
573  return dataset_offset;
574 }
575 
576 
577 
578 template <int dim>
580  const Quadrature<dim - 1> &quadrature)
581 {
582  return project_to_all_faces(ReferenceCells::get_hypercube<dim>(), quadrature);
583 }
584 
585 
586 template <int dim>
589  const Quadrature<dim - 1> &quadrature)
590 {
591  return project_to_all_faces(reference_cell,
592  hp::QCollection<dim - 1>(quadrature));
593 }
594 
595 
596 /* -------------- declaration of explicit specializations ------------- */
597 
598 #ifndef DOXYGEN
599 
600 
601 template <>
602 void
604  const unsigned int,
605  std::vector<Point<1>> &);
606 template <>
607 void
609  const Quadrature<0> &,
610  const unsigned int,
611  std::vector<Point<1>> &);
612 template <>
613 void
615  const unsigned int face_no,
616  std::vector<Point<2>> &q_points);
617 template <>
618 void
619 QProjector<2>::project_to_face(const ReferenceCell reference_cell,
620  const Quadrature<1> & quadrature,
621  const unsigned int face_no,
622  std::vector<Point<2>> &q_points);
623 template <>
624 void
626  const unsigned int face_no,
627  std::vector<Point<3>> &q_points);
628 template <>
629 void
630 QProjector<3>::project_to_face(const ReferenceCell reference_cell,
631  const Quadrature<2> & quadrature,
632  const unsigned int face_no,
633  std::vector<Point<3>> &q_points);
634 
635 template <>
638  const hp::QCollection<0> &quadrature);
639 
640 
641 template <>
642 void
644  const unsigned int,
645  const unsigned int,
646  std::vector<Point<1>> &,
647  const RefinementCase<0> &);
648 template <>
649 void
651  const Quadrature<0> &,
652  const unsigned int,
653  const unsigned int,
654  std::vector<Point<1>> &,
655  const RefinementCase<0> &);
656 template <>
657 void
659  const unsigned int face_no,
660  const unsigned int subface_no,
661  std::vector<Point<2>> &q_points,
662  const RefinementCase<1> &);
663 template <>
664 void
666  const Quadrature<1> & quadrature,
667  const unsigned int face_no,
668  const unsigned int subface_no,
669  std::vector<Point<2>> &q_points,
670  const RefinementCase<1> &);
671 template <>
672 void
674  const unsigned int face_no,
675  const unsigned int subface_no,
676  std::vector<Point<3>> & q_points,
677  const RefinementCase<2> &face_ref_case);
678 template <>
679 void
681  const Quadrature<2> & quadrature,
682  const unsigned int face_no,
683  const unsigned int subface_no,
684  std::vector<Point<3>> & q_points,
685  const RefinementCase<2> &face_ref_case);
686 
687 template <>
690 template <>
693  const Quadrature<0> &quadrature);
694 
695 
696 #endif // DOXYGEN
698 
699 #endif
static const unsigned int invalid_unsigned_int
Definition: types.h:196
static Quadrature< dim > project_to_child(const Quadrature< dim > &quadrature, const unsigned int child_no)
Definition: qprojector.cc:1238
const unsigned int dataset_offset
Definition: qprojector.h:530
static void project_to_subface(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_faces(const Quadrature< dim - 1 > &quadrature)
Definition: qprojector.h:579
static DataSetDescriptor cell()
Definition: qprojector.h:563
static Quadrature< dim > project_to_all_subfaces(const SubQuadrature &quadrature)
static Quadrature< dim > project_to_line(const Quadrature< 1 > &quadrature, const Point< dim > &p1, const Point< dim > &p2)
Definition: qprojector.cc:1323
void reference_cell(const ReferenceCell &reference_cell, Triangulation< dim, spacedim > &tria)
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:394
static void project_to_face(const SubQuadrature &quadrature, const unsigned int face_no, std::vector< Point< dim >> &q_points)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:393
static DataSetDescriptor subface(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)
#define DEAL_II_DEPRECATED
Definition: config.h:156
static Quadrature< dim > project_to_all_children(const Quadrature< dim > &quadrature)
Definition: qprojector.cc:1282
static DataSetDescriptor face(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:1365