Reference documentation for deal.II version 9.4.1
\(\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\}}\)
Loading...
Searching...
No Matches
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
33
85template <int dim>
86class Quadrature : public Subscriptor
87{
88public:
94 using SubQuadrature = Quadrature<dim == 0 ? 0 : dim - 1>;
95
104 explicit Quadrature(const unsigned int n_quadrature_points = 0);
105
114 Quadrature(const SubQuadrature &, const Quadrature<1> &);
115
135 explicit Quadrature(const Quadrature<dim != 1 ? 1 : 0> &quadrature_1d);
136
140 Quadrature(const Quadrature<dim> &q);
141
146 Quadrature(Quadrature<dim> &&) noexcept = default;
147
153 Quadrature(const std::vector<Point<dim>> &points,
154 const std::vector<double> & weights);
155
163 Quadrature(const std::vector<Point<dim>> &points);
164
169 Quadrature(const Point<dim> &point);
170
174 virtual ~Quadrature() override = default;
175
180 Quadrature &
181 operator=(const Quadrature<dim> &);
182
187 Quadrature &
188 operator=(Quadrature<dim> &&) = default; // NOLINT
189
193 bool
194 operator==(const Quadrature<dim> &p) const;
195
200 void
201 initialize(const std::vector<Point<dim>> &points,
202 const std::vector<double> & weights);
203
207 unsigned int
208 size() const;
209
213 const Point<dim> &
214 point(const unsigned int i) const;
215
219 const std::vector<Point<dim>> &
220 get_points() const;
221
225 double
226 weight(const unsigned int i) const;
227
231 const std::vector<double> &
232 get_weights() const;
233
238 std::size_t
239 memory_consumption() const;
240
246 template <class Archive>
247 void
248 serialize(Archive &ar, const unsigned int version);
249
255 bool
257
276#ifndef DOXYGEN
277 typename std::conditional<dim == 1,
278 std::array<Quadrature<1>, dim>,
279 const std::array<Quadrature<1>, dim> &>::type
280#else
281 const std::array<Quadrature<1>, dim> &
282#endif
283 get_tensor_basis() const;
284
285protected:
290 std::vector<Point<dim>> quadrature_points;
291
296 std::vector<double> weights;
297
306
311 std::unique_ptr<std::array<Quadrature<1>, dim>> tensor_basis;
312};
313
314
323template <int dim>
324class QAnisotropic : public Quadrature<dim>
325{
326public:
331 QAnisotropic(const Quadrature<1> &qx);
332
336 QAnisotropic(const Quadrature<1> &qx, const Quadrature<1> &qy);
337
341 QAnisotropic(const Quadrature<1> &qx,
342 const Quadrature<1> &qy,
343 const Quadrature<1> &qz);
344};
345
346
370template <int dim>
371class QIterated : public Quadrature<dim>
372{
373public:
380 QIterated(const Quadrature<1> &base_quadrature, const unsigned int n_copies);
381
393 QIterated(const Quadrature<1> & base_quadrature,
394 const std::vector<Point<1>> &intervals);
395
400 "The quadrature formula you provided cannot be used "
401 "as the basis for iteration.");
402};
403
404
405
408#ifndef DOXYGEN
409
410// ------------------- inline and template functions ----------------
411
412
413template <int dim>
414inline unsigned int
416{
417 return weights.size();
418}
419
420
421template <int dim>
422inline const Point<dim> &
423Quadrature<dim>::point(const unsigned int i) const
424{
425 AssertIndexRange(i, size());
426 return quadrature_points[i];
427}
428
429
430
431template <int dim>
432double
433Quadrature<dim>::weight(const unsigned int i) const
434{
435 AssertIndexRange(i, size());
436 return weights[i];
437}
438
439
440
441template <int dim>
442inline const std::vector<Point<dim>> &
444{
445 return quadrature_points;
446}
447
448
449
450template <int dim>
451inline const std::vector<double> &
453{
454 return weights;
455}
456
457
458
459template <int dim>
460inline bool
462{
463 return is_tensor_product_flag;
464}
465
466
467
468template <int dim>
469template <class Archive>
470inline void
471Quadrature<dim>::serialize(Archive &ar, const unsigned int)
472{
473 // forward to serialization
474 // function in the base class.
475 ar &static_cast<Subscriptor &>(*this);
476
477 ar &quadrature_points &weights;
478}
479
480
481
482/* -------------- declaration of explicit specializations ------------- */
483
484template <>
485Quadrature<0>::Quadrature(const unsigned int);
486template <>
488 const Quadrature<1> &);
489template <>
491template <>
493
494template <>
496
497template <>
499
500template <>
501QIterated<1>::QIterated(const Quadrature<1> &base_quadrature,
502 const unsigned int n_copies);
503
504#endif // DOXYGEN
506
507#endif
Definition: point.h:111
QIterated(const Quadrature< 1 > &base_quadrature, const unsigned int n_copies)
Definition: quadrature.cc:630
std::vector< Point< dim > > quadrature_points
Definition: quadrature.h:290
std::unique_ptr< std::array< Quadrature< 1 >, dim > > tensor_basis
Definition: quadrature.h:311
bool is_tensor_product() const
std::size_t memory_consumption() const
Definition: quadrature.cc:313
const Point< dim > & point(const unsigned int i) const
bool is_tensor_product_flag
Definition: quadrature.h:305
Quadrature(const unsigned int n_quadrature_points=0)
Definition: quadrature.cc:42
const std::array< Quadrature< 1 >, dim > & get_tensor_basis() const
Definition: quadrature.cc:325
double weight(const unsigned int i) const
const std::vector< double > & get_weights() const
void initialize(const std::vector< Point< dim > > &points, const std::vector< double > &weights)
Definition: quadrature.cc:52
Quadrature(Quadrature< dim > &&) noexcept=default
std::vector< double > weights
Definition: quadrature.h:296
const std::vector< Point< dim > > & get_points() const
unsigned int size() const
void serialize(Archive &ar, const unsigned int version)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
static ::ExceptionBase & ExcInvalidQuadratureFormula()
#define AssertIndexRange(index, range)
Definition: exceptions.h:1732
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:487
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
STL namespace.