Reference documentation for deal.II version GIT f6a5d312c9 2023-10-04 08:50: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 - 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_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 
197  Quadrature(std::vector<Point<dim>> &&points, std::vector<double> &&weights);
198 
206  Quadrature(const std::vector<Point<dim>> &points);
207 
212  Quadrature(const Point<dim> &point);
213 
217  virtual ~Quadrature() override = default;
218 
223  Quadrature &
224  operator=(const Quadrature<dim> &);
225 
230  Quadrature &
231  operator=(Quadrature<dim> &&) = default; // NOLINT
232 
236  bool
237  operator==(const Quadrature<dim> &p) const;
238 
243  void
244  initialize(const std::vector<Point<dim>> &points,
245  const std::vector<double> &weights);
246 
250  unsigned int
251  size() const;
252 
256  bool
257  empty() const;
258 
262  const Point<dim> &
263  point(const unsigned int i) const;
264 
268  const std::vector<Point<dim>> &
269  get_points() const;
270 
274  double
275  weight(const unsigned int i) const;
276 
280  const std::vector<double> &
281  get_weights() const;
282 
287  std::size_t
288  memory_consumption() const;
289 
295  template <class Archive>
296  void
297  serialize(Archive &ar, const unsigned int version);
298 
304  bool
306 
325 #ifndef DOXYGEN
326  typename std::conditional<dim == 1,
327  std::array<Quadrature<1>, dim>,
328  const std::array<Quadrature<1>, dim> &>::type
329 #else
330  const std::array<Quadrature<1>, dim> &
331 #endif
332  get_tensor_basis() const;
333 
334 protected:
339  std::vector<Point<dim>> quadrature_points;
340 
345  std::vector<double> weights;
346 
355 
360  std::unique_ptr<std::array<Quadrature<1>, dim>> tensor_basis;
361 };
362 
363 
372 template <int dim>
373 class QAnisotropic : public Quadrature<dim>
374 {
375 public:
380  QAnisotropic(const Quadrature<1> &qx);
381 
385  QAnisotropic(const Quadrature<1> &qx, const Quadrature<1> &qy);
386 
390  QAnisotropic(const Quadrature<1> &qx,
391  const Quadrature<1> &qy,
392  const Quadrature<1> &qz);
393 };
394 
395 
419 template <int dim>
420 class QIterated : public Quadrature<dim>
421 {
422 public:
429  QIterated(const Quadrature<1> &base_quadrature, const unsigned int n_copies);
430 
442  QIterated(const Quadrature<1> &base_quadrature,
443  const std::vector<Point<1>> &intervals);
444 
449  "The quadrature formula you provided cannot be used "
450  "as the basis for iteration.");
451 };
452 
453 
454 
457 #ifndef DOXYGEN
458 
459 // ------------------- inline and template functions ----------------
460 
461 
462 template <int dim>
463 inline unsigned int
464 Quadrature<dim>::size() const
465 {
466  return weights.size();
467 }
468 
469 
470 
471 template <int dim>
472 inline bool
474 {
475  return weights.empty();
476 }
477 
478 
479 
480 template <int dim>
481 inline const Point<dim> &
482 Quadrature<dim>::point(const unsigned int i) const
483 {
484  AssertIndexRange(i, size());
485  return quadrature_points[i];
486 }
487 
488 
489 
490 template <int dim>
491 double
492 Quadrature<dim>::weight(const unsigned int i) const
493 {
494  AssertIndexRange(i, size());
495  return weights[i];
496 }
497 
498 
499 
500 template <int dim>
501 inline const std::vector<Point<dim>> &
503 {
504  return quadrature_points;
505 }
506 
507 
508 
509 template <int dim>
510 inline const std::vector<double> &
512 {
513  return weights;
514 }
515 
516 
517 
518 template <int dim>
519 inline bool
521 {
522  return is_tensor_product_flag;
523 }
524 
525 
526 
527 template <int dim>
528 template <class Archive>
529 inline void
530 Quadrature<dim>::serialize(Archive &ar, const unsigned int)
531 {
532  // forward to serialization
533  // function in the base class.
534  ar &static_cast<Subscriptor &>(*this);
535 
536  ar &quadrature_points &weights;
537 }
538 
539 
540 
541 /* -------------- declaration of explicit specializations ------------- */
542 
543 template <>
544 Quadrature<0>::Quadrature(const unsigned int);
545 template <>
547  const Quadrature<1> &);
548 template <>
550 template <>
552 
553 template <>
555 
556 template <>
558 
559 template <>
560 QIterated<1>::QIterated(const Quadrature<1> &base_quadrature,
561  const unsigned int n_copies);
562 
563 #endif // DOXYGEN
565 
566 #endif
Definition: point.h:112
QAnisotropic(const Quadrature< 1 > &qx)
Definition: quadrature.cc:366
QIterated(const Quadrature< 1 > &base_quadrature, const unsigned int n_copies)
Definition: quadrature.cc:639
std::vector< Point< dim > > quadrature_points
Definition: quadrature.h:339
std::unique_ptr< std::array< Quadrature< 1 >, dim > > tensor_basis
Definition: quadrature.h:360
const std::vector< Point< dim > > & get_points() const
bool is_tensor_product() const
std::size_t memory_consumption() const
Definition: quadrature.cc:326
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:354
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:338
const Point< dim > & point(const unsigned int i) const
Quadrature(Quadrature< dim > &&) noexcept=default
std::vector< double > weights
Definition: quadrature.h:345
bool empty() const
unsigned int size() const
void serialize(Archive &ar, const unsigned int version)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
#define AssertIndexRange(index, range)
Definition: exceptions.h:1857
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:490
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:490