Reference documentation for deal.II version GIT d2cc530c04 2023-03-22 15:10: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\}}\)
quadrature.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2022 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_quadrature_h
17 #define dealii_quadrature_h
18 
19 
20 #include <deal.II/base/config.h>
21 
22 #include <deal.II/base/point.h>
24 
25 #include <array>
26 #include <memory>
27 #include <vector>
28 
30 
121 template <int dim>
122 class Quadrature : public Subscriptor
123 {
124 public:
130  using SubQuadrature = Quadrature<dim == 0 ? 0 : dim - 1>;
131 
140  explicit Quadrature(const unsigned int n_quadrature_points = 0);
141 
150  Quadrature(const SubQuadrature &, const Quadrature<1> &);
151 
171  explicit Quadrature(const Quadrature<dim != 1 ? 1 : 0> &quadrature_1d);
172 
176  Quadrature(const Quadrature<dim> &q);
177 
182  Quadrature(Quadrature<dim> &&) noexcept = default;
183 
189  Quadrature(const std::vector<Point<dim>> &points,
190  const std::vector<double> & weights);
191 
199  Quadrature(const std::vector<Point<dim>> &points);
200 
205  Quadrature(const Point<dim> &point);
206 
210  virtual ~Quadrature() override = default;
211 
216  Quadrature &
217  operator=(const Quadrature<dim> &);
218 
223  Quadrature &
224  operator=(Quadrature<dim> &&) = default; // NOLINT
225 
229  bool
230  operator==(const Quadrature<dim> &p) const;
231 
236  void
237  initialize(const std::vector<Point<dim>> &points,
238  const std::vector<double> & weights);
239 
243  unsigned int
244  size() const;
245 
249  const Point<dim> &
250  point(const unsigned int i) const;
251 
255  const std::vector<Point<dim>> &
256  get_points() const;
257 
261  double
262  weight(const unsigned int i) const;
263 
267  const std::vector<double> &
268  get_weights() const;
269 
274  std::size_t
275  memory_consumption() const;
276 
282  template <class Archive>
283  void
284  serialize(Archive &ar, const unsigned int version);
285 
291  bool
293 
312 #ifndef DOXYGEN
313  typename std::conditional<dim == 1,
314  std::array<Quadrature<1>, dim>,
315  const std::array<Quadrature<1>, dim> &>::type
316 #else
317  const std::array<Quadrature<1>, dim> &
318 #endif
319  get_tensor_basis() const;
320 
321 protected:
326  std::vector<Point<dim>> quadrature_points;
327 
332  std::vector<double> weights;
333 
342 
347  std::unique_ptr<std::array<Quadrature<1>, dim>> tensor_basis;
348 };
349 
350 
359 template <int dim>
360 class QAnisotropic : public Quadrature<dim>
361 {
362 public:
367  QAnisotropic(const Quadrature<1> &qx);
368 
372  QAnisotropic(const Quadrature<1> &qx, const Quadrature<1> &qy);
373 
377  QAnisotropic(const Quadrature<1> &qx,
378  const Quadrature<1> &qy,
379  const Quadrature<1> &qz);
380 };
381 
382 
406 template <int dim>
407 class QIterated : public Quadrature<dim>
408 {
409 public:
416  QIterated(const Quadrature<1> &base_quadrature, const unsigned int n_copies);
417 
429  QIterated(const Quadrature<1> & base_quadrature,
430  const std::vector<Point<1>> &intervals);
431 
436  "The quadrature formula you provided cannot be used "
437  "as the basis for iteration.");
438 };
439 
440 
441 
444 #ifndef DOXYGEN
445 
446 // ------------------- inline and template functions ----------------
447 
448 
449 template <int dim>
450 inline unsigned int
451 Quadrature<dim>::size() const
452 {
453  return weights.size();
454 }
455 
456 
457 template <int dim>
458 inline const Point<dim> &
459 Quadrature<dim>::point(const unsigned int i) const
460 {
461  AssertIndexRange(i, size());
462  return quadrature_points[i];
463 }
464 
465 
466 
467 template <int dim>
468 double
469 Quadrature<dim>::weight(const unsigned int i) const
470 {
471  AssertIndexRange(i, size());
472  return weights[i];
473 }
474 
475 
476 
477 template <int dim>
478 inline const std::vector<Point<dim>> &
480 {
481  return quadrature_points;
482 }
483 
484 
485 
486 template <int dim>
487 inline const std::vector<double> &
489 {
490  return weights;
491 }
492 
493 
494 
495 template <int dim>
496 inline bool
498 {
499  return is_tensor_product_flag;
500 }
501 
502 
503 
504 template <int dim>
505 template <class Archive>
506 inline void
507 Quadrature<dim>::serialize(Archive &ar, const unsigned int)
508 {
509  // forward to serialization
510  // function in the base class.
511  ar &static_cast<Subscriptor &>(*this);
512 
513  ar &quadrature_points &weights;
514 }
515 
516 
517 
518 /* -------------- declaration of explicit specializations ------------- */
519 
520 template <>
521 Quadrature<0>::Quadrature(const unsigned int);
522 template <>
524  const Quadrature<1> &);
525 template <>
527 template <>
529 
530 template <>
532 
533 template <>
535 
536 template <>
537 QIterated<1>::QIterated(const Quadrature<1> &base_quadrature,
538  const unsigned int n_copies);
539 
540 #endif // DOXYGEN
542 
543 #endif
Definition: point.h:110
QAnisotropic(const Quadrature< 1 > &qx)
Definition: quadrature.cc:353
QIterated(const Quadrature< 1 > &base_quadrature, const unsigned int n_copies)
Definition: quadrature.cc:626
std::vector< Point< dim > > quadrature_points
Definition: quadrature.h:326
std::unique_ptr< std::array< Quadrature< 1 >, dim > > tensor_basis
Definition: quadrature.h:347
const std::vector< Point< dim > > & get_points() const
bool is_tensor_product() const
std::size_t memory_consumption() const
Definition: quadrature.cc:313
void initialize(const std::vector< Point< dim >> &points, const std::vector< double > &weights)
Definition: quadrature.cc:52
bool is_tensor_product_flag
Definition: quadrature.h:341
Quadrature(const unsigned int n_quadrature_points=0)
Definition: quadrature.cc:42
const std::vector< double > & get_weights() const
double weight(const unsigned int i) const
const std::array< Quadrature< 1 >, dim > & get_tensor_basis() const
Definition: quadrature.cc:325
const Point< dim > & point(const unsigned int i) const
Quadrature(Quadrature< dim > &&) noexcept=default
std::vector< double > weights
Definition: quadrature.h:332
unsigned int size() const
void serialize(Archive &ar, const unsigned int version)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:474
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:475
#define AssertIndexRange(index, range)
Definition: exceptions.h:1827
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:488
static ::ExceptionBase & ExcInvalidQuadratureFormula()
void quadrature_points(const Triangulation< dim, spacedim > &triangulation, const Quadrature< dim > &quadrature, const std::vector< std::vector< BoundingBox< spacedim >>> &global_bounding_boxes, ParticleHandler< dim, spacedim > &particle_handler, const Mapping< dim, spacedim > &mapping=(ReferenceCells::get_hypercube< dim >() .template get_default_linear_mapping< dim, spacedim >()), const std::vector< std::vector< double >> &properties={})
Definition: generators.cc:493