Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
qprojector.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2005 - 2018 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/geometry_info.h>
21 #include <deal.II/base/quadrature.h>
22 
23 DEAL_II_NAMESPACE_OPEN
24 
25 
28 
29 
75 template <int dim>
77 {
78 public:
83  using SubQuadrature = Quadrature<dim - 1>;
84 
90  static void
91  project_to_face(const SubQuadrature & quadrature,
92  const unsigned int face_no,
93  std::vector<Point<dim>> &q_points);
94 
100  static Quadrature<dim>
101  project_to_face(const SubQuadrature &quadrature, const unsigned int face_no);
102 
112  static void
113  project_to_subface(const SubQuadrature & quadrature,
114  const unsigned int face_no,
115  const unsigned int subface_no,
116  std::vector<Point<dim>> & q_points,
117  const RefinementCase<dim - 1> &ref_case =
119 
129  static Quadrature<dim>
130  project_to_subface(const SubQuadrature & quadrature,
131  const unsigned int face_no,
132  const unsigned int subface_no,
133  const RefinementCase<dim - 1> &ref_case =
135 
152  static Quadrature<dim>
153  project_to_all_faces(const SubQuadrature &quadrature);
154 
168  static Quadrature<dim>
169  project_to_all_subfaces(const SubQuadrature &quadrature);
170 
181  static Quadrature<dim>
182  project_to_child(const Quadrature<dim> &quadrature,
183  const unsigned int child_no);
184 
194  static Quadrature<dim>
195  project_to_all_children(const Quadrature<dim> &quadrature);
196 
201  static Quadrature<dim>
202  project_to_line(const Quadrature<1> &quadrature,
203  const Point<dim> & p1,
204  const Point<dim> & p2);
205 
219  {
220  public:
227 
234  static DataSetDescriptor
235  cell();
236 
247  static DataSetDescriptor
248  face(const unsigned int face_no,
249  const bool face_orientation,
250  const bool face_flip,
251  const bool face_rotation,
252  const unsigned int n_quadrature_points);
253 
266  static DataSetDescriptor
267  subface(const unsigned int face_no,
268  const unsigned int subface_no,
269  const bool face_orientation,
270  const bool face_flip,
271  const bool face_rotation,
272  const unsigned int n_quadrature_points,
273  const internal::SubfaceCase<dim> ref_case =
275 
282  operator unsigned int() const;
283 
284  private:
288  const unsigned int dataset_offset;
289 
294  DataSetDescriptor(const unsigned int dataset_offset);
295  };
296 
297 private:
305  static Quadrature<2>
306  reflect(const Quadrature<2> &q);
307 
317  static Quadrature<2>
318  rotate(const Quadrature<2> &q, const unsigned int n_times);
319 };
320 
324 // ------------------- inline and template functions ----------------
325 
326 
327 
328 template <int dim>
330  const unsigned int dataset_offset)
331  : dataset_offset(dataset_offset)
332 {}
333 
334 
335 template <int dim>
337  : dataset_offset(numbers::invalid_unsigned_int)
338 {}
339 
340 
341 
342 template <int dim>
345 {
346  return 0;
347 }
348 
349 
350 
351 template <int dim>
353 {
354  return dataset_offset;
355 }
356 
357 
358 /* -------------- declaration of explicit specializations ------------- */
359 
360 #ifndef DOXYGEN
361 
362 
363 template <>
364 void
366  const unsigned int,
367  std::vector<Point<1>> &);
368 template <>
369 void
371  const unsigned int face_no,
372  std::vector<Point<2>> &q_points);
373 template <>
374 void
376  const unsigned int face_no,
377  std::vector<Point<3>> &q_points);
378 
379 template <>
382 
383 
384 template <>
385 void
387  const unsigned int,
388  const unsigned int,
389  std::vector<Point<1>> &,
390  const RefinementCase<0> &);
391 template <>
392 void
394  const unsigned int face_no,
395  const unsigned int subface_no,
396  std::vector<Point<2>> &q_points,
397  const RefinementCase<1> &);
398 template <>
399 void
401  const unsigned int face_no,
402  const unsigned int subface_no,
403  std::vector<Point<3>> & q_points,
404  const RefinementCase<2> &face_ref_case);
405 
406 template <>
409 
410 template <>
411 QIterated<1>::QIterated(const Quadrature<1> &base_quadrature,
412  const unsigned int n_copies);
413 
414 
415 #endif // DOXYGEN
416 DEAL_II_NAMESPACE_CLOSE
417 
418 #endif
static Quadrature< dim > project_to_child(const Quadrature< dim > &quadrature, const unsigned int child_no)
Definition: quadrature.cc:1063
const unsigned int dataset_offset
Definition: qprojector.h:288
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< 2 > rotate(const Quadrature< 2 > &q, const unsigned int n_times)
Definition: quadrature.cc:450
static DataSetDescriptor cell()
Definition: qprojector.h:344
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: quadrature.cc:1116
static Quadrature< 2 > reflect(const Quadrature< 2 > &q)
Definition: quadrature.cc:432
static Quadrature< dim > project_to_all_faces(const SubQuadrature &quadrature)
static void project_to_face(const SubQuadrature &quadrature, const unsigned int face_no, std::vector< Point< dim >> &q_points)
QIterated(const Quadrature< 1 > &base_quadrature, const unsigned int n_copies)
Definition: quadrature.cc:1777
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)
static Quadrature< dim > project_to_all_children(const Quadrature< dim > &quadrature)
Definition: quadrature.cc:1090
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: quadrature.cc:1139