Reference documentation for deal.II version GIT db2fd67796 2022-12-07 23:35: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 
87 template <int dim>
88 class Quadrature : public Subscriptor
89 {
90 public:
96  using SubQuadrature = Quadrature<dim == 0 ? 0 : dim - 1>;
97 
106  explicit Quadrature(const unsigned int n_quadrature_points = 0);
107 
116  Quadrature(const SubQuadrature &, const Quadrature<1> &);
117 
137  explicit Quadrature(const Quadrature<dim != 1 ? 1 : 0> &quadrature_1d);
138 
142  Quadrature(const Quadrature<dim> &q);
143 
148  Quadrature(Quadrature<dim> &&) noexcept = default;
149 
155  Quadrature(const std::vector<Point<dim>> &points,
156  const std::vector<double> & weights);
157 
165  Quadrature(const std::vector<Point<dim>> &points);
166 
171  Quadrature(const Point<dim> &point);
172 
176  virtual ~Quadrature() override = default;
177 
182  Quadrature &
183  operator=(const Quadrature<dim> &);
184 
189  Quadrature &
190  operator=(Quadrature<dim> &&) = default; // NOLINT
191 
195  bool
196  operator==(const Quadrature<dim> &p) const;
197 
202  void
203  initialize(const std::vector<Point<dim>> &points,
204  const std::vector<double> & weights);
205 
209  unsigned int
210  size() const;
211 
215  const Point<dim> &
216  point(const unsigned int i) const;
217 
221  const std::vector<Point<dim>> &
222  get_points() const;
223 
227  double
228  weight(const unsigned int i) const;
229 
233  const std::vector<double> &
234  get_weights() const;
235 
240  std::size_t
241  memory_consumption() const;
242 
248  template <class Archive>
249  void
250  serialize(Archive &ar, const unsigned int version);
251 
257  bool
259 
278 #ifndef DOXYGEN
279  typename std::conditional<dim == 1,
280  std::array<Quadrature<1>, dim>,
281  const std::array<Quadrature<1>, dim> &>::type
282 #else
283  const std::array<Quadrature<1>, dim> &
284 #endif
285  get_tensor_basis() const;
286 
287 protected:
292  std::vector<Point<dim>> quadrature_points;
293 
298  std::vector<double> weights;
299 
308 
313  std::unique_ptr<std::array<Quadrature<1>, dim>> tensor_basis;
314 };
315 
316 
325 template <int dim>
326 class QAnisotropic : public Quadrature<dim>
327 {
328 public:
333  QAnisotropic(const Quadrature<1> &qx);
334 
338  QAnisotropic(const Quadrature<1> &qx, const Quadrature<1> &qy);
339 
343  QAnisotropic(const Quadrature<1> &qx,
344  const Quadrature<1> &qy,
345  const Quadrature<1> &qz);
346 };
347 
348 
372 template <int dim>
373 class QIterated : public Quadrature<dim>
374 {
375 public:
382  QIterated(const Quadrature<1> &base_quadrature, const unsigned int n_copies);
383 
395  QIterated(const Quadrature<1> & base_quadrature,
396  const std::vector<Point<1>> &intervals);
397 
402  "The quadrature formula you provided cannot be used "
403  "as the basis for iteration.");
404 };
405 
406 
407 
410 #ifndef DOXYGEN
411 
412 // ------------------- inline and template functions ----------------
413 
414 
415 template <int dim>
416 inline unsigned int
417 Quadrature<dim>::size() const
418 {
419  return weights.size();
420 }
421 
422 
423 template <int dim>
424 inline const Point<dim> &
425 Quadrature<dim>::point(const unsigned int i) const
426 {
427  AssertIndexRange(i, size());
428  return quadrature_points[i];
429 }
430 
431 
432 
433 template <int dim>
434 double
435 Quadrature<dim>::weight(const unsigned int i) const
436 {
437  AssertIndexRange(i, size());
438  return weights[i];
439 }
440 
441 
442 
443 template <int dim>
444 inline const std::vector<Point<dim>> &
446 {
447  return quadrature_points;
448 }
449 
450 
451 
452 template <int dim>
453 inline const std::vector<double> &
455 {
456  return weights;
457 }
458 
459 
460 
461 template <int dim>
462 inline bool
464 {
465  return is_tensor_product_flag;
466 }
467 
468 
469 
470 template <int dim>
471 template <class Archive>
472 inline void
473 Quadrature<dim>::serialize(Archive &ar, const unsigned int)
474 {
475  // forward to serialization
476  // function in the base class.
477  ar &static_cast<Subscriptor &>(*this);
478 
479  ar &quadrature_points &weights;
480 }
481 
482 
483 
484 /* -------------- declaration of explicit specializations ------------- */
485 
486 template <>
487 Quadrature<0>::Quadrature(const unsigned int);
488 template <>
490  const Quadrature<1> &);
491 template <>
493 template <>
495 
496 template <>
498 
499 template <>
501 
502 template <>
503 QIterated<1>::QIterated(const Quadrature<1> &base_quadrature,
504  const unsigned int n_copies);
505 
506 #endif // DOXYGEN
508 
509 #endif
Definition: point.h:111
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:292
std::unique_ptr< std::array< Quadrature< 1 >, dim > > tensor_basis
Definition: quadrature.h:313
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:307
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:298
unsigned int size() const
void serialize(Archive &ar, const unsigned int version)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:458
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:459
#define AssertIndexRange(index, range)
Definition: exceptions.h:1760
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:487
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