BarycentricRational.hpp
1 // Distributed under the MIT License.
2 // See LICENSE.txt for details.
3 
4 #pragma once
5 
6 #include <cstddef>
7 #include <unistd.h> // IWYU pragma: keep
8 #include <vector>
9 
10 namespace PUP {
11 class er;
12 } // namespace PUP
13 
14 namespace intrp {
15 /*!
16  * \ingroup NumericalAlgorithmsGroup
17  * \brief A barycentric rational interpolation class
18  *
19  * The class builds a barycentric rational interpolant of a specified order
20  * using the `x_values` and `y_values` passed into the constructor.
21  * Barycentric interpolation requires \f$3N\f$ storage, and costs
22  * \f$\mathcal{O}(Nd)\f$ to construct, where \f$N\f$ is the size of the x- and
23  * y-value vectors and \f$d\f$ is the order of the interpolant. The evaluation
24  * cost is \f$\mathcal{O}(N)\f$ compared to \f$\mathcal{O}(d)\f$ of a spline
25  * method, but constructing the barycentric interpolant does not require any
26  * derivatives of the function to be known.
27  *
28  * The interpolation function is
29  *
30  * \f[
31  * \mathcal{I}(x)=\frac{\sum_{i=0}^{N-1}w_i y_i /
32  * (x-x_i)}{\sum_{i=0}^{N-1}w_i/(x-x_i)}
33  * \f]
34  *
35  * where \f$w_i\f$ are the weights. The weights are computed using
36  *
37  * \f[
38  * w_k=\sum_{i=k-d\\0\le i < N-d}^{k}(-1)^{i}
39  * \prod_{j=i\\j\ne k}^{i+d}\frac{1}{x_k-x_j} \f]
40  *
41  * \requires `x_values.size() == y_values.size()` and
42  * `x_values_.size() >= order`
43  */
45  public:
46  BarycentricRational() noexcept = default;
48  std::vector<double> y_values, size_t order) noexcept;
49 
50  double operator()(double x_to_interp_to) const noexcept;
51 
52  size_t order() const noexcept;
53 
54  // clang-tidy: no runtime references
55  void pup(PUP::er& p) noexcept; // NOLINT
56 
57  private:
58  std::vector<double> x_values_;
59  std::vector<double> y_values_;
60  std::vector<double> weights_;
61  ssize_t order_{0};
62 };
63 } // namespace intrp
Definition: Strahlkorper.hpp:14
Definition: AddTemporalIdsToInterpolationTarget.hpp:17
A barycentric rational interpolation class.
Definition: BarycentricRational.hpp:44