Interpolations.jl 中的 MATLAB 样式 griddedInterpolant 用法

Question

我正在尝试将我的 MATLAB 程序转换为 Julia。该程序的一个关键特性是使用 MATLAB 中的 griddedInterpolant 函数。我找到了 Julia 替代品 (Interpolations.jl)，并且我在二维案例中做了一个简单的测试，以确保我理解它是如何工作的。不过，这个特定的程序使用 4-D 数组，我似乎无法弄清楚如何在 2 维之外使用 Interpolations.jl 方法。

一个非常简化的 MATLAB 行为示例。请注意，这是一个简单的示例，因为 V 是平坦的，但您明白了。我的程序要做的是在循环中更改V、kp1gv 和kp2gv 的值。

nK_1 = 50;
nK_2 = 50;
nKP_1 = 10;
nKP_2 = 10;
K1 = linspace(0,100,nK_1);
K2 = linspace(0,100,nK_2);
KP1 = linspace(0,100,nKP_1);
KP2 = linspace(0,100,nKP_1);

nX = 3;
nY = 3;
X = [0.9,1,1.1];
Y = [0.95,1,1.05];

[k1gv,k2gv,xgv,ygv] = ndgrid(K1,K2,X,Y);  % creates grid vectors
V = ones(nK_1,nK_2,nX,nY);    % value func to be interpolated
Fit = griddedInterpolant(k1gv,k2gv,xgv,ygv,V,'linear');   % fitted val fun


[kp1gv,kp2gv,xpgv,ypgv,kk1,kk2] = ndgrid(KP1,KP2,X,Y,K1,K2);
Fitted = Fit(kp1gv,kp2gv,xpgv,ypgv);

我正在寻找的是在 Julia 中复制它，或者完全像在 MATLAB 中一样，或者通过重写它来更快。这是我现在的尝试：

nK_1 = 50
nK_2 = 50
nKP_1 = 10
nKP_2 = 10
K1 = linspace(0,100,nK_1)
K2 = linspace(0,100,nK_2)
KP1 = linspace(0,100,nKP_1)
KP2 = linspace(0,100,nKP_1)

nX = 3
nY = 3
X = [0.9,1,1.1]
Y = [0.95,1,1.05]

V = ones(nK_1,nK_2,nX,nY)

# using the more compact notation here, would be similar to
# griddedInterpolant({K1,K2,X,Y},V,'linear') -- no issues with this fit
Fit = interpolate((K1,K2,X,Y),V,Gridded(Linear()))

# this is how I'm getting ndgrid functionality
KP1gv = Float64[i for i in KP1, j in KP2, x in X, y in Y, a in K1, b in K2]
KP2gv = Float64[j for i in KP1, j in KP2, x in X, y in Y, a in K1, b in K2]
Xgv = Float64[x for i in KP1, j in KP2, x in X, y in Y, a in K1, b in K2]
Ygv = Float64[y for i in KP1, j in KP2, x in X, y in Y, a in K1, b in K2]

# what I want to work
Fitted = Fit[KP1gv,KP2gv,Xgv,Ygv]

我已阅读有关 Interpolations.jl 的文档，我知道这应该是可行的，但我似乎无法让它工作。

Answer 1

与其构造所有ndgrid 坐标数组，不如这样做：

fitted = [Fit[i,j,k,l] for i in KP1, j in KP2, k in X, l in Y]

但有几个提示：

• 任何对性能至关重要的东西（如上面的行）都应该放在一个函数中。这绝对是至关重要的。请参阅手册的performance tips 页面。

• 这似乎没有很好的记录，但Gridded 是为您有一个矩形网格但间距可能不规则的情况而设计的。在规则间距的情况下，您可以做得更好

  itp = interpolate(V, Linear(), OnGrid()) sitp = scale(itp, K1, K2, x, y)

其中x 和y 也是linspace 对象。