34 std::unique_ptr<Function>
splineInterpolation(
const std::vector<double>& x,
const std::vector<double>& y,
41 std::vector<double> h (n, 0.);
42 for(
int i=0; i<n; i++)
45 std::vector<double> mu (n, 0.);
46 std::vector<double> z (n+1, 0.);
48 for (
int i=1; i<n; ++i) {
49 g = 2.*(x[i+1]-x[i-1]) - h[i-1]*mu[i-1];
51 z[i] = (3.*(y[i+1]*h[i-1] - y[i]*(x[i+1]-x[i-1])+ y[i-1]*h[i]) / (h[i-1]*h[i]) - h[i-1] * z[i-1]) /
g;
55 std::vector<double> a (n, 0.);
56 std::vector<double> b (n, 0.);
57 std::vector<double> c (n+1, 0.);
58 std::vector<double> d (n, 0.);
63 for (
int j=n-1; j>=0; j--) {
65 c[j] = z[j] - mu[j] * c[j+1];
66 b[j] = (y[j+1] - y[j]) / h[j] - h[j] * (c[j+1] + 2. * c[j]) / 3.;
67 d[j] = (c[j+1] - c[j]) / (3. * h[j]);
72 for (
int i=0; i<n; i++) {
74 double x_2 = x_1 * x_1;
75 double x_3 = x_1 * x_2;
76 a[i] = a[i] + b[i]*x_1 + c[i]*x_2 + d[i]*x_3;
77 b[i] = b[i] + 2.*c[i]*x_1 + 3.*d[i]*x_2;
78 c[i] = c[i] + 3.*d[i]*x_1;
82 std::vector<std::unique_ptr<Function>> functions;
84 for (
int i=0; i<n; i++) {
85 functions.emplace_back(std::unique_ptr<Function>(
new Polynomial{{a[i],b[i],c[i],d[i]}}));
89 std::vector<double> x_copy(x);
90 x_copy.front() = std::numeric_limits<double>::lowest();
91 x_copy.back() = std::numeric_limits<double>::max();
92 return std::unique_ptr<Function>(
new Piecewise{x_copy, std::move(functions)});
95 return std::unique_ptr<Function>(
new Piecewise{x, std::move(functions)});
Represents a polynomial function.
std::unique_ptr< Function > splineInterpolation(const std::vector< double > &x, const std::vector< double > &y, bool extrapolate)
Performs cubic spline interpolation for the given set of data points.
Represents a piecewise function.