Применение динамической функции к каждой точке в фазовом пространстве (представленной 2D-матрицей)

Question

У меня есть матрица целых чисел, phase_space формы (n,n), где каждая запись представляет количество точек в этом месте в пространстве. У меня также есть две матрицы обновления u_x, u_y также формы (n,n), с целыми числами в диапазоне 0,...,n, указывающими, где моя динамическая система занимает каждую соответствующую точку в пространстве. Я хочу "применить" матрицы обновления к фазовому пространству итеративно.

Например, если

>>>u_x
array([[1, 2, 1],
       [0, 1, 2],
       [0, 0, 0]])
>>>u_y
array([[2, 1, 2],
       [1, 0, 1],
       [2, 2, 0]])
>>>phase_space 
array([[1, 1, 1],
       [1, 1, 1],
       [1, 1, 1]])

Я хочу

>>>new_phase_space
array([[1., 1., 2.],
       [1., 0., 2.],
       [0., 2., 0.]])

Мой текущий (рабочий)Решение состоит в следующем:

for i in range(n):
    for j in range(n):
        new_phase_space[u_x[i, j], u_y[i, j]] += phase_space[i,j] 

Есть ли способ векторизовать это?

Answer 1

Вы можете использовать pandas.DataFrame.groupby() для накопления всех ходов с одинаковыми координатами в phase_space:

new_phase_space + (pd.DataFrame(phase_space)
           .stack()
           .groupby([u_x.ravel(), u_y.ravel()])
           .sum()
           .unstack(fill_value=0)
           .values
)

Выход:

array([[2., 2., 4.],
       [2., 0., 4.],
       [0., 4., 0.]])

Answer 2

Мы можем использовать np.bincount -

M,N = u_x.max()+1,u_y.max()+1
ids = u_x*N+u_y
out = np.bincount(ids.ravel(),phase_space.ravel(),minlength=M*N).reshape(M,N)

Образец запуска на более общей установке -

In [14]: u_x
Out[14]: 
array([[1, 2, 1],
       [0, 1, 4],
       [0, 0, 0]])

In [15]: u_y
Out[15]: 
array([[2, 1, 2],
       [6, 0, 1],
       [2, 6, 0]])

In [17]: phase_space
Out[17]: 
array([[1, 1, 1],
       [5, 1, 1],
       [1, 1, 1]])

In [18]: out
Out[18]: 
array([[1., 0., 1., 0., 0., 0., 6.],
       [1., 0., 2., 0., 0., 0., 0.],
       [0., 1., 0., 0., 0., 0., 0.],
       [0., 0., 0., 0., 0., 0., 0.],
       [0., 1., 0., 0., 0., 0., 0.]])

Мы также могли бы использовать разреженные матрицыособенно если речь идет о памяти -

from scipy.sparse import csr_matrix,coo_matrix

out = coo_matrix( (phase_space.ravel(), (u_x.ravel(), u_y.ravel())), shape = (M,N))

Вывод будет разреженной матрицей. Чтобы преобразовать в плотный, используйте out.toarray().