如何使用Python Scipy正确实现3D矢量场插值

核心错误点 #

• 低级变量误用：插值函数中对v、w分量reshape时错误传入了u变量，导致三个矢量分量全部和u分量一致，完全丢失了原始v、w的场特征
• 网格维度顺序不匹配：numpy.meshgrid默认使用indexing="xy"，生成的原始场数组维度顺序为(y轴长度, x轴长度, z轴长度)，但你reshape时使用的顺序是(x长度, y长度, z长度)，同时RegularGridInterpolator要求坐标元组顺序和数据数组的轴顺序严格对应，二者不匹配导致插值结果完全错位

修复后的完整插值函数 #

def interpolate_field(x,y,z,u,v,w,new_points):
    x_unique = np.unique(x)
    y_unique = np.unique(y)
    z_unique = np.unique(z)
    # 按meshgrid生成时的维度顺序reshape，原始indexing=xy对应shape (len(y), len(x), len(z))
    u = np.reshape(u, (len(y_unique), len(x_unique), len(z_unique)))
    v = np.reshape(v, (len(y_unique), len(x_unique), len(z_unique)))
    w = np.reshape(w, (len(y_unique), len(x_unique), len(z_unique)))
    # 坐标元组顺序和数据轴顺序严格对应：第一轴对应y，第二对应x，第三对应z
    u_int_f = RegularGridInterpolator((y_unique, x_unique, z_unique), u)
    v_int_f = RegularGridInterpolator((y_unique, x_unique, z_unique), v)
    w_int_f = RegularGridInterpolator((y_unique, x_unique, z_unique), w)
    
    u_int = u_int_f(new_points)
    v_int = v_int_f(new_points)
    w_int = w_int_f(new_points)
    
    return u_int, v_int, w_int

可选优化方案 #

如果不想调整维度顺序，可以在原始生成网格时指定indexing="ij"，生成的数组维度顺序就会和(x长度, y长度, z长度)完全匹配，后续reshape和插值时不需要调整坐标顺序，代码逻辑更直观。

内容的提问来源于stack exchange，提问作者henry