Anonymous user
Thiele's interpolation formula: Difference between revisions
Updated D entry
(Updated D entry) |
|||
Line 355:
struct Domain {
const real b, e, s;
auto range() const /*pure nothrow*/ {
}
}
Line 365 ⟶ 364:
real eval0(alias RY, alias X, alias Y)(in real x) pure nothrow {
real a = 0.0L;
foreach_reverse (immutable i; 2 .. X.length - 3)
a = (x - X[i]) / (RY[i] - RY[i-2] + a);
return Y[1] + (x - X[1]) / (RY[1] + a);
}
immutable struct Thiele {
immutable real[] Y, X, rhoY, rhoX;
this(real[] y, real[] x)
in {
assert(x.length > 2, "at leat 3 values");
assert(x.length == y.length, "input arrays not of same size");
} body {
this.Y = y.idup; // Not nothrow.
this.X = x.idup;
rhoY = rhoN(Y, X);
Line 385 ⟶ 384:
this(in real delegate(real) f,
Domain d=Domain(0.0L, 1.55L, 0.05L))
immutable /*pure nothrow*/ {
auto xrng = d.range().array();▼
}
auto rhoN(immutable(real)[] y, immutable(real)[] x)
pure /*nothrow*/ {
auto p = new real[
p[
p[1][0 .. $ - 1]
foreach (immutable int j; 2 .. N - 1) {
immutable M = N - j - 1;
p[j][0..M] = p[j-2][1..M+1] + (x[0..M] - x[j..M+j]) /
(p[j-1][0 .. M] - p[j-1][1 .. M+1]);
}
return p.map!q{ a[1] }
}
alias eval = eval0!(rhoY, X, Y)
alias inverse = eval0!(rhoX, Y, X)
}
void main() {
// can't pass sin,cos,tan directly
writefln(" %d interpolating points\n", tsin.X.length);
writefln("std.math.sin(0.5): %20.18f",
writefln(" Thiele sin(0.5): %20.18f\n", tsin.eval(0.5L));
|