我正在慢慢地从C转向Python。这次我需要从给定的网格中计算出偏导数。我知道如何在C语言中执行此操作,因此目前我仅使用内联适配器,即
def dz(x,X,Y,Z,dx):
y = numpy.zeros((X,Y,Z), dtype='double');
code = """
int i, j, k;
for (i=0; i<X-1; i++){
for(k=0; k<Y; k++){
for (j=0; j<Z; j++){
y[i,k,j] = (x[i+1, k, j] - x[i, k, j])/dx;
}
}
}
for (j=0; j<Z; j++){
for(k=0; k<Y; k++){
y[X-1,k,j] = - x[X-1, k, j]/dx;
}
}
"""
weave.inline(code, ['x', 'y', 'dx', 'X', 'Y', 'Z'], \
type_converters=converters.blitz, compiler = 'gcc');
return y;
其中x
和y
是3D numpy数组,如您所见,第二个循环代表边界条件。当然,我可以在纯Python中实现相同的逻辑,但是代码效率很低。但是,我想知道是否可以使用纯numpy计算偏导数?任何人都可以提供任何帮助,我将不胜感激。
np.diff
可能是最惯用的numpy方法:
np.diff
您可能也对y = np.empty_like(x)
y[:-1] = np.diff(x, axis=0) / dx
y[-1] = -x[-1] / dx
感兴趣,尽管此函数对输入数组的所有维度(而不是单个维度)进行渐变。
y = np.empty_like(x)
y[:-1] = (x[1:] - x[:-1]) / dx
y[-1] = -x[-1] / dx
y = np.empty_like(x)
y[:, :-1] = (x[:, 1:] - x[:, :-1]) / dx
y[:, -1] = -x[:, -1] / dx
现在也支持沿单个方向求导数。
def dz(x,dx):
y = numpy.zeros(x.shape, dtype='double')
y[:-1] = (x[1:] - x[:-1]) / dx
y[-1] = -x[-1] / dx
return y
将给出沿轴numpy.gradient
的单个偏导数:
numpy.gradient
参数a = np.array([[1,2,3],[2,3,5]])
np.gradient(a, axis=0)
指定一组求导数的方向。