• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
1[/
2Copyright (c) 2020 Nick Thompson
3Use, modification and distribution are subject to the
4Boost Software License, Version 1.0. (See accompanying file
5LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
6]
7
8[section:makima Modified Akima interpolation]
9
10[heading Synopsis]
11
12    #include <boost/math/interpolators/makima.hpp>
13
14    namespace boost::math::interpolators {
15
16    template <class RandomAccessContainer>
17    class makima
18    {
19    public:
20
21        using Real = RandomAccessContainer::value_type;
22
23        makima(RandomAccessContainer&& abscissas, RandomAccessContainer&& ordinates,
24               Real left_endpoint_derivative = std::numeric_limits<Real>::quiet_NaN(),
25               Real right_endpoint_derivative = std::numeric_limits<Real>::quiet_NaN());
26
27        Real operator()(Real x) const;
28
29        Real prime(Real x) const;
30
31        void push_back(Real x, Real y);
32
33        friend std::ostream& operator<<(std::ostream & os, const makima & m);
34    };
35
36    } // namespaces
37
38
39[heading Modified Akima Interpolation]
40
41The modified Akima interpolant takes non-equispaced data and interpolates between them via cubic Hermite polynomials whose slopes are chosen by a modification of a geometric construction proposed by [@https://doi.org/10.1145/321607.321609 Akima].
42The modification is given by [@https://blogs.mathworks.com/cleve/2019/04/29/makima-piecewise-cubic-interpolation/ Cosmin Ionita] and agrees with Matlab's version.
43The interpolant is /C/[super 1] and evaluation has [bigo](log(/N/)) complexity.
44This is faster than barycentric rational interpolation, but also less smooth.
45An example usage is as follows:
46
47    std::vector<double> x{1, 5, 9 , 12};
48    std::vector<double> y{8,17, 4, -3};
49    using boost::math::interpolators::makima;
50    auto spline = makima(std::move(x), std::move(y));
51    // evaluate at a point:
52    double z = spline(3.4);
53    // evaluate derivative at a point:
54    double zprime = spline.prime(3.4);
55
56Periodically, it is helpful to see what data the interpolator has, and the slopes it has chosen.
57This can be achieved via
58
59    std::cout << spline << "\n";
60
61Note that the interpolator is pimpl'd, so that copying the class is cheap, and hence it can be shared between threads.
62(The call operator and `.prime()` are threadsafe.)
63
64One unique aspect of this interpolator is that it can be updated in constant time.
65Hence we can use `boost::circular_buffer` to do real-time interpolation:
66
67    #include <boost/circular_buffer.hpp>
68    ...
69    boost::circular_buffer<double> initial_x{1,2,3,4};
70    boost::circular_buffer<double> initial_y{4,5,6,7};
71    auto circular_akima = makima(std::move(initial_x), std::move(initial_y));
72    // interpolate via call operation:
73    double y = circular_akima(3.5);
74    // add new data:
75    circular_akima.push_back(5, 8);
76    // interpolat at 4.5:
77    y = circular_akima(4.5);
78
79
80
81[$../graphs/makima_vs_cubic_b.svg]
82
83The modified Akima spline compared to the cubic /B/-spline.
84The modified Akima spline oscillates less than the cubic spline, but has less smoothness and is not exact on quadratic polynomials.
85
86
87[heading Complexity and Performance]
88
89The complexity and performance is identical to that of the cubic Hermite interpolator, since this object simply constructs derivatives and forwards the data to `cubic_hermite.hpp`.
90
91[endsect]
92[/section:makima]
93