以矩阵作为输出在 python 中的多维数组中进行插值

Question

在我的问题中，我有一个根据一个或多个参数变化并输出矩阵的函数。现在对该函数进行采样。所以假设函数 F 依赖于 x 和 y，那么对于一些 x-y 网格，我为每个网格点都有一个矩阵。

在 MATLAB 中，我使用 4D 矩阵解决此问题，其中前两个维度是矩阵本身，其余维度是 x-y 网格的坐标。请参阅以下示例：

n = 2;
% x-y grid 
xg = [-1, 1, 7];
sz_xg = size(xg,2);
yg = [3, 4, 9];
sz_yg = size(yg,2);

% some function:
fcn = @(x,y) [3*x+sin(y), cos(x); exp(-x), 5+y]; % some n x n matrix function

% Define grid
A = zeros(n,n,sz_xg,sz_yg);
for i = 1:sz_xg
    for j = 1:sz_yg
        A(1:n, 1:n, i, j) = fcn(xg(i), yg(j));
    end
end

xpoint = 5;
ypoint = 6;

Ainterp = interpn(1:n, 1:n, xg, yg, A, 1:n, 1:n, xpoint, ypoint)

我也尝试在 Python 中修复此问题，但没有成功。当我尝试实现完全相同的示例时，我得到一个 2x1 向量...

我的尝试：

import numpy as np
from scipy.interpolate import interpn

n = 2
xg = np.array([-1, 1, 7])
yg = np.array([3, 4, 9])

# some function
def fcn(x,y):
    return np.array([[3*x+np.sin(y), np.cos(x)],[np.exp(-x), 5+y]])

A = np.zeros((n,n,xg.size,yg.size))
for i in range(xg.size):
    for j in range(yg.size):
        A[:,:,i,j]=fcn(xg[i],yg[j])

xpoint = 5
ypoint = 6

At = interpn((np.arange(n), np.arange(n), xg, yg),
             A,(np.arange(0,n),np.arange(0,n),xpoint,ypoint))

你能帮帮我吗？

Answer 1

interpn 的实现在 Matlab 和 scipy 之间略有不同。

Matlab returns 输入向量的排列，因此当您为第一维和第二维指定 [1,2] 时，它 returns 的值在

[(1,1) (1,2);
 (2,1) (2,2)]

但是scipy interpn不是这种情况，scipy使用并行数组来确定插值的位置，这意味着第一个输出是从两个数组中的第一个输入构造的，并且第二个输出由两个数组中的第二个输入构成，依此类推。

为了获得您想要的行为，您必须使用下面描述的更有效的方法传递上面列出的整个矩阵：

grid_points = np.meshgrid(np.arange(0,n),np.arange(0,n),sparse=True)

At = interpn((np.arange(n), np.arange(n), xg, yg),
             A,(grid_points[0],grid_points[1],xpoint,ypoint))

基本上，meshgrid(sparse=True) 将以高效和内存压缩的方式为您生成我上面写的矩阵。

在我看来，python 方法比 matlab 更可预测，因为它更容易确定输出数组的维数，这与 matlab 不同，后者的输出维数取决于输入不是单个值，这使得很难预测输出中维度的顺序。

OP 注释： 当插值的输出是向量时，您可能会出错。因此，例如，At 应该是 3x1，并且它被插入 3 个变量，网格为 5。因此，A 将是 3x1x5x5x5。那么在误差计算中很可能会出现错误。我通过在 interpn 函数

之前添加来避免这种情况

np.seterr(divide='ignore', invalid='ignore')

和interpn函数之后

np.seterr(divide='warn', invalid='warn')

以矩阵作为输出在 python 中的多维数组中进行插值

Interpolate in multidimensional array in python with matrix as output

python

matlab

interpolation

matrix