Так что, в основном, я пытаюсь сделать конечно-разностное на 2d массиве, не делая слишком много для циклов. Я хотел бы иметь гессенскую матрицу массива и градиент. Поэтому мне нужны как производные первого и второго порядка массива.
Это может быть достигнуто путем вычисления следующего уравнения в массиве. 
Чтобы справиться с границами, мы вычисляем его только для внутренних точек, поэтому код для этого производного может выглядеть примерно так:
arr = np.random.rand(16).reshape(4,4)
result = np.zeros_like(arr)
w, h = arr.shape
for i in range(1, w-1):
for j in range(1, h-1):
result[i,j] = (arr[i+1, j] - arr[i-1, j]) / (2*dx)
Это дает правильный ответ, но он может быть очень медленным по сравнению с nu numpy операциями, поэтому я подумал про себя. По сути, это просто свертка с ядром, которое выглядит следующим образом
kernel = [1, 0 , -1]
Поэтому мы выполняем следующий код
from scipy.sigmal import convolve
result = np.pad((convolve(arr,kernel,mode='same',
method = 'direct')/(2*dx))[1:-1, 1:-1], 1).T
Поскольку мы имеем дело только с внутренними точками, мы их обрезаем и дополнить нулями впоследствии, чтобы подражать тому, что произошло в предыдущем наивном случае.
Это работает! Но с некоторыми массивами, среднеквадратическая ошибка между наивным случаем и свёртыми небесными ракетами. Таким образом, кажется, что числовая ошибка очень сильно увеличивается в некоторых случаях.
Мне бы хотелось, чтобы скорость, полученная в результате свертки, была стабильной из наивного случая. Любая помощь?