mir.interpolate.linear

Linear!(F, N, X) linear(F, size_t N = 1, X = F)(Repeat!(N, Slice!(RCI!(immutable(X)))) grid, Slice!(RCI!(const(F)), N) values);

Constructs multivariate linear interpolant with nodes on rectilinear grid.

Parameters:

Repeat!(N, Slice!(RCI!(immutable(X)))) grid	x values for interpolant
Slice!(RCI!(const(F)), N) values	f(x) values for interpolant

Constraints grid, values must have the same length >= 2

Returns:

Linear

Examples:

R -> R: Linear interpolation

import mir.algorithm.iteration;
import mir.ndslice;
import mir.math.common: approxEqual;

static immutable x = [0, 1, 2, 3, 5.00274, 7.00274, 10.0055, 20.0137, 30.0192];
static immutable y = [0.0011, 0.0011, 0.0030, 0.0064, 0.0144, 0.0207, 0.0261, 0.0329, 0.0356,];
static immutable xs = [1, 2, 3, 4.00274, 5.00274, 6.00274, 7.00274, 8.00548, 9.00548, 10.0055, 11.0055, 12.0082, 13.0082, 14.0082, 15.0082, 16.011, 17.011, 18.011, 19.011, 20.0137, 21.0137, 22.0137, 23.0137, 24.0164, 25.0164, 26.0164, 27.0164, 28.0192, 29.0192, 30.0192];

auto interpolant = linear!double(x.rcslice!(immutable double), y.rcslice!(const double));

static immutable data = [0.0011, 0.0030, 0.0064, 0.0104, 0.0144, 0.0176, 0.0207, 0.0225, 0.0243, 0.0261, 0.0268, 0.0274, 0.0281, 0.0288, 0.0295, 0.0302, 0.0309, 0.0316, 0.0322, 0.0329, 0.0332, 0.0335, 0.0337, 0.0340, 0.0342, 0.0345, 0.0348, 0.0350, 0.0353, 0.0356];

assert(xs.sliced.vmap(interpolant).all!((a, b) => approxEqual(a, b, 1e-4, 1e-4))(data));

Examples:

R^2 -> R: Bilinear interpolation

import mir.math.common: approxEqual;
import mir.ndslice;
alias appreq = (a, b) => approxEqual(a, b, 10e-10, 10e-10);

//// set test function ////
enum y_x0 = 2;
enum y_x1 = -7;
enum y_x0x1 = 3;

// this function should be approximated very well
alias f = (x0, x1) => y_x0 * x0 + y_x1 * x1 + y_x0x1 * x0 * x1 - 11;

///// set interpolant ////
static immutable x0 = [-1.0, 2, 8, 15];
static immutable x1 = [-4.0, 2, 5, 10, 13];

auto grid = cartesian(x0, x1)
    .map!f
    .rcslice
    .lightConst;

auto interpolant = 
    linear!(double, 2)(
        x0.rcslice!(immutable double),
        x1.rcslice!(immutable double),
        grid
    );

///// compute test data ////
auto test_grid = cartesian(x0.sliced + 1.23, x1.sliced + 3.23);
auto real_data = test_grid.map!f;
auto interp_data = test_grid.vmap(interpolant);
///// verify result ////
assert(all!appreq(interp_data, real_data));

//// check derivatives ////
auto z0 = 1.23;
auto z1 = 3.21;
auto d = interpolant.withDerivative(z0, z1);
assert(appreq(interpolant(z0, z1), f(z0, z1)));
assert(appreq(d[0][0], f(z0, z1)));
assert(appreq(d[1][0], y_x0 + y_x0x1 * z1));
assert(appreq(d[0][1], y_x1 + y_x0x1 * z0));
assert(appreq(d[1][1], y_x0x1));

Examples:

R^3 -> R: Trilinear interpolation

import mir.math.common: approxEqual;
import mir.ndslice;
alias appreq = (a, b) => approxEqual(a, b, 10e-10, 10e-10);

///// set test function ////
enum y_x0 = 2;
enum y_x1 = -7;
enum y_x2 = 5;
enum y_x0x1 = 10;
enum y_x0x1x2 = 3;

// this function should be approximated very well
static auto f(double x0, double x1, double x2)
{
    return y_x0 * x0 + y_x1 * x1 + y_x2 * x2 + y_x0x1 * x0 * x1 + y_x0x1x2 * x0 * x1 * x2 - 11;
}

///// set interpolant ////
static immutable x0 = [-1.0, 2, 8, 15];
static immutable x1 = [-4.0, 2, 5, 10, 13];
static immutable x2 = [3, 3.7, 5];
auto grid = cartesian(x0, x1, x2)
    .map!f
    .as!(const double)
    .rcslice;

auto interpolant = linear!(double, 3)(
        x0.rcslice!(immutable double),
        x1.rcslice!(immutable double),
        x2.rcslice!(immutable double),
        grid);

///// compute test data ////
auto test_grid = cartesian(x0.sliced + 1.23, x1.sliced + 3.23, x2.sliced - 3);
auto real_data = test_grid.map!f;
auto interp_data = test_grid.vmap(interpolant);
///// verify result ////
assert(all!appreq(interp_data, real_data));

//// check derivatives ////
auto z0 = 1.23;
auto z1 = 3.21;
auto z2 = 4;
auto d = interpolant.withDerivative(z0, z1, z2);
assert(appreq(interpolant(z0, z1, z2), f(z0, z1, z2)));
assert(appreq(d[0][0][0], f(z0, z1, z2)));
assert(appreq(d[1][0][0], y_x0 + y_x0x1 * z1 + y_x0x1x2 * z1 * z2));
assert(appreq(d[0][1][0], y_x1 + y_x0x1 * z0 + y_x0x1x2 * z0 * z2));
assert(appreq(d[1][1][0], y_x0x1 + y_x0x1x2 * z2));
assert(appreq(d[1][1][1], y_x0x1x2));

struct Linear(F, size_t N = 1, X = F) if (N && (N <= 6));

Multivariate linear interpolant with nodes on rectilinear grid.

Slice!(RCI!(const(F)), N) _data;

Aligned buffer allocated with mir.internal.memory. For internal use.

Repeat!(N, RCI!(immutable(X))) _grid;

Grid iterators. For internal use.

@nogc @safe this(Repeat!(N, Slice!(RCI!(immutable(X)))) grid, Slice!(RCI!(const(F)), N) data);

const @property Linear lightConst()();

const @property scope Slice!(RCI!(immutable(X))) grid(size_t dimension = 0)() return
if (dimension < N);

const @property scope @trusted immutable(X)[] gridScopeView(size_t dimension = 0)() return
if (dimension < N);

const @property scope size_t intervalCount(size_t dimension = 0)();

Returns:

intervals count.

const @property scope size_t[N] gridShape()();

enum uint derivativeOrder;

template opCall(uint derivative = 0) if (derivative <= derivativeOrder)

const scope @trusted auto opCall(X...)(in X xs)
if (X.length == N);

(x) operator.

Complexity O(log(grid.length))

alias withDerivative = opCall!1;

struct LinearKernel(X);

const @property auto lightConst()();

immutable @property auto lightImmutable()();

this()(X x0, X x1, X x);

template opCall(uint derivative = 0) if (derivative <= 1)

auto opCall(Y)(in Y y0, in Y y1);

alias withDerivative = opCall!1;