Multidimensional linear interpolation

New to Julia
,
1

I have some 2D Array data x that I want to (linearly) interpolate along a single vector a uniformly in each dimension (with the Interpolations package). An example is below:

using Interpolations
a   = collect(range(0, stop=5, length=10))
x   = vcat([collect(range(0, stop=5i, length=10)) for i in 1:3]'...)

but neither of the two below attempts works

itp = interpolate((a,a,a), x, Gridded(Linear()))
itp = interpolate((a,), x, Gridded(Linear()))

I could create an array of interpolants for each dimension manually

itp = [interpolate((a,), x[i,:], Gridded(Linear())) for i in 1:3]

but was wondering if I could efficiently do this in “one go.”

2

Please include what packages you have usinged in your MWE.

We don’t know what Gridded(Linear()) you mean, its not in the standard library.
So it must be from some package, but there are over 2000 julia packages.

:slight_smile:

4 Likes
3

It seems to me like you don’t actually want multivariate linear interpolation. Bivariate gridded interpolation is for when you have two vectors x and y and a matrixZ=F(x,y) for each point the combined x,y space (these are the “knots”).

itp =interpolate((x,y),Z,Gridded(Linear()))

It sounds to me like you want to take some irregularly spaced data and put it on a grid a. In this case you do need to interpolate each time independently (using the actual knots) and then evaluate the interpolant at a. Interpolations.jl is pretty quick though, especially for univariate linear interpolation, so you shouldn’t worry about it being too slow.

Did I misunderstand?

2 Likes
4

Thanks @tbeason and @oxinabox, I added more detail (using Interpolations). Basically I have a function f: \mathbb{R} \to \mathbb{R}^3 evaluated at [f(a_1), f(a_2), \ldots, f(a_N)] and I want to interpolate linearly f(a) between the a_i.

5

f : \mathbb{R} \rightarrow \mathbb{R}^3 … can’t say I’ve used interpolations in a context like that before. It sounds to me like you’d need to do it separately.

6

I’m a little confused about what you were trying above, but if your goal is

then it should work out of the box if you do exactly what the math tells you to do:

using Interpolations, StaticArrays
f(a) = SVector{3}(2a, 3a, 4a)  # replace with your real function
agrid = 0:0.5:5
vals = [f(a) for a in agrid]
itp = interpolate((agrid,), vals, Gridded(Linear()))

julia> itp(0.25)
3-element SArray{Tuple{3},Float64,1,3} with indices SOneTo(3):
 0.5 
 0.75
 1.0

julia> f(0.25)    # since the function happens to be linear this should match to within numeric precision
3-element SArray{Tuple{3},Float64,1,3} with indices SOneTo(3):
 0.5 
 0.75
 1.0

The confusion might stem from representing vals as a 2d array rather than a 1d array.

If your grid is evenly spaced, performance-wise you’ll be much better off using a scaled BSpline than Gridded.

4 Likes
7

Ah I see, I wasn’t clear, sorry. Instead of a function f, I have evaluations of f on a grid, and I’d like to interpolate linearly on the grid without knowing / having access to f itself.

8

f was just for demonstration purposes. If you have agrid and vals from another source the itp = interpolate(...) call is the same regardless.

9

Ah, the difference is that I have vals as a 2D array.

10

Right. Given the math you posted above, you will likely be better off representing it as an array-of-arrays. (In particular the fact that it’s 3-dimensional is encoded statically, which often allows for some efficiencies.) If you must store it as a 2d array, try itp = interpolate((1:3, agrid,), data, (NoInterp(), Gridded(Linear()))).

2 Likes
11

Cool answer. I also have a similar question with the question owner then I fortunately find your solution.
I have a question that why StaticArrays can help Interpolations on interpolating vector value function? How does it work?