Reference documentation for deal.II version 9.5.0
\(\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 - 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
121template <int dim>
123{
124public:
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 const Point<dim> &
257 point(const unsigned int i) const;
258
262 const std::vector<Point<dim>> &
263 get_points() const;
264
268 double
269 weight(const unsigned int i) const;
270
274 const std::vector<double> &
275 get_weights() const;
276
281 std::size_t
282 memory_consumption() const;
283
289 template <class Archive>
290 void
291 serialize(Archive &ar, const unsigned int version);
292
298 bool
300
319#ifndef DOXYGEN
320 typename std::conditional<dim == 1,
321 std::array<Quadrature<1>, dim>,
322 const std::array<Quadrature<1>, dim> &>::type
323#else
324 const std::array<Quadrature<1>, dim> &
325#endif
326 get_tensor_basis() const;
327
328protected:
333 std::vector<Point<dim>> quadrature_points;
334
339 std::vector<double> weights;
340
349
354 std::unique_ptr<std::array<Quadrature<1>, dim>> tensor_basis;
355};
356
357
366template <int dim>
367class QAnisotropic : public Quadrature<dim>
368{
369public:
374 QAnisotropic(const Quadrature<1> &qx);
375
379 QAnisotropic(const Quadrature<1> &qx, const Quadrature<1> &qy);
380
384 QAnisotropic(const Quadrature<1> &qx,
385 const Quadrature<1> &qy,
386 const Quadrature<1> &qz);
387};
388
389
413template <int dim>
414class QIterated : public Quadrature<dim>
415{
416public:
423 QIterated(const Quadrature<1> &base_quadrature, const unsigned int n_copies);
424
436 QIterated(const Quadrature<1> & base_quadrature,
437 const std::vector<Point<1>> &intervals);
438
443 "The quadrature formula you provided cannot be used "
444 "as the basis for iteration.");
445};
446
447
448
451#ifndef DOXYGEN
452
453// ------------------- inline and template functions ----------------
454
455
456template <int dim>
457inline unsigned int
459{
460 return weights.size();
461}
462
463
464template <int dim>
465inline const Point<dim> &
466Quadrature<dim>::point(const unsigned int i) const
467{
468 AssertIndexRange(i, size());
469 return quadrature_points[i];
470}
471
472
473
474template <int dim>
475double
476Quadrature<dim>::weight(const unsigned int i) const
477{
478 AssertIndexRange(i, size());
479 return weights[i];
480}
481
482
483
484template <int dim>
485inline const std::vector<Point<dim>> &
487{
488 return quadrature_points;
489}
490
491
492
493template <int dim>
494inline const std::vector<double> &
496{
497 return weights;
498}
499
500
501
502template <int dim>
503inline bool
505{
506 return is_tensor_product_flag;
507}
508
509
510
511template <int dim>
512template <class Archive>
513inline void
514Quadrature<dim>::serialize(Archive &ar, const unsigned int)
515{
516 // forward to serialization
517 // function in the base class.
518 ar &static_cast<Subscriptor &>(*this);
519
520 ar &quadrature_points &weights;
521}
522
523
524
525/* -------------- declaration of explicit specializations ------------- */
526
527template <>
528Quadrature<0>::Quadrature(const unsigned int);
529template <>
531 const Quadrature<1> &);
532template <>
534template <>
536
537template <>
539
540template <>
542
543template <>
544QIterated<1>::QIterated(const Quadrature<1> &base_quadrature,
545 const unsigned int n_copies);
546
547#endif // DOXYGEN
549
550#endif
Definition point.h:112
QIterated(const Quadrature< 1 > &base_quadrature, const unsigned int n_copies)
std::vector< Point< dim > > quadrature_points
Definition quadrature.h:333
std::unique_ptr< std::array< Quadrature< 1 >, dim > > tensor_basis
Definition quadrature.h:354
bool is_tensor_product() const
std::size_t memory_consumption() const
const Point< dim > & point(const unsigned int i) const
bool is_tensor_product_flag
Definition quadrature.h:348
Quadrature(const unsigned int n_quadrature_points=0)
Definition quadrature.cc:42
const std::array< Quadrature< 1 >, dim > & get_tensor_basis() const
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:339
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:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
static ::ExceptionBase & ExcInvalidQuadratureFormula()
#define AssertIndexRange(index, range)
#define DeclExceptionMsg(Exception, defaulttext)
Definition exceptions.h:488
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={})
STL namespace.