Non-uniform interpolation

[TheBB]

I'm opening this issue based on this discussion. I've looked at your source code relatively thoroughly, and I think I understand how it works, even though it's a bit more meta than I'm used to.

Am I correct in assuming that non-uniform interpolation is best represented as a new subtype of GridRepresentation, for example with knot vector as argument to the constructor?

If so, I'd like to start working on that, beginning with constants and working my way up.

[tomasaschan]

Wonderful! =)

I'm not exactly sure exactly how to best tackle non-uniform grids, actually. Depending on how different the math behind the interpolations is, it might be better to either build on the current infrastructure (and add new subtypes to the current configuration flags) or build something new alongside it.

I'm thinking that the API has to be similar to the one for uniform interpolation, so where we currently do

itp = interpolate(yvalues, configs...) # configs is all the configuration type objects

with an assumed x-grid of 1:size(yvalues,dim) in each dimension, maybe we want something like

itp = interpolate(xgrids, yvalues, configs...)

for non-uniform interpolation. xgrids here is tricky; we want it to be general enough to work in any dimensionality, without constraints on the number of grid points in different directions (other than that they should match those of yvalues, of course...). To be space-efficient, it would be preferable to specify the grid in any one dimension as a Vector, so maybe a tuple of Vectors or something similar is appropriate.

Since much of the current infrastructure is focused on solving the interpolation problem in 1D and then generalizing through recursion, maybe the same approach is good enough to start with for non-uniforms as well, in which case xgrids and yvalues are simply Vectors, both of them, and we can figure out API details later.

[timholy]

@TheBB, are your non-uniform grids cartesian product grids, or could you have knots at random points? If the former, I like the tuple-of-Vectors interface, with a special case for a single-vector input (for 1d).

I assume you already understand this (you seem to catch on quick!), but just in case: by encoding something about the order of your interpolation into the type parameters, you'll be able to generate efficient code that, for example, doesn't require for loops.

[timholy]

Also, if we want to target julia 0.4, keep in mind that stagedfunctions are available (if they would be helpful).

[tomasaschan]

Re. Julia version: I have seen 0.4 as the main target for this library from the start, so feel free to use stuff that requires it. I won't be testing on 0.4 until 0.4-pre is released (mainly because I'm lazy...) but that shouldn't stop anyone else =)

[TheBB]

@TheBB, are your non-uniform grids cartesian product grids, or could you have knots at random points?

Goodness, no. Just tensor product grids.

[timholy]

Whew 😄.

[NateNesbitt]

Thanks for your work on this package - I'm excited for the non-uniform interpolation to be available! I'm new to GitHub, is this the proper place for a comment like this, or is there somewhere more appropriate for expressing interest in features like the package you're working on here?

Best,

Nate

[TheBB]

Here is fine. Unfortunately I had to start working on something rather different the last couple of months, so I haven't been able to dedicate time to this. :-(

[pwl]

I don't know if this is of any help but Fornberg [1] provides an algorithm for approximating derivatives (of any order) of a function given on a non-uniform 1d mesh. In particular, it can be used to compute a value of a function (which is the zeroth-order derivative) at any given point (so not only at a mesh point). The algorithm is pretty simple to implement and understand/modify, although I don't know if his approach can be generalized to tensor product grids.

[1] http://www.jstor.org/stable/2653239

[tomasaschan]

@pwl Thanks, that's a good reference. I don't have access to the article atm (the full article isn't available for free) but I can probably get access through my university. I will probably focus my efforts on getting the B-spline functionality (see #5) working, and stabilize the API (see #31), but I'll keep this in mind for the future if no-one else gets there before me :)

[tomasaschan]

I'm wondering if this can be considered solved in terms of Gridded (#56)...