Я использую np.einsum
для расчета потока материала в графике (от 1 узла до 4 узлов в этом примере).Количество потока задается как amount
(amount.shape == (1, 1, 2)
размеры определяют определенные критерии, назовем их a
, b
, c
).
Булева матрица route
определяетдопустимый поток на основе критериев a
, b
, c
в y
(route.shape == (4, 1, 1, 2)
; yabc
).Я помечаю размеры y
, a
, b
, c
.abc
эквивалентны amount
s измерениям abc
, y
- направление потока (0, 1, 2 или 3).Чтобы определить количество материала в y
, я вычисляю np.einsum('abc,yabc->y', amount, route)
и получаю вектор y-dim с потоками в y
.Есть также неявная расстановка приоритетов маршрута.Например, любой route[0, ...] == True
равен False
для любого y=1..3
, любой route[1, ...] == True
равен False
для следующих более высоких маршрутов y-dim и так далее.route[3, ...]
(последний индекс y) определяет маршрут перехвата, то есть его значения равны True
, когда предыдущими значениями индекса y были Ложь ((route[0] ^ route[1] ^ route[2] ^ route[3]).all() == True
).
Это прекрасно работает.Однако, когда я ввожу другой критерий (измерение) x
, который существует только в route
, но не в amount
, эта логика, кажется, нарушается.Приведенный ниже код демонстрирует проблему:
>>> import numpy as np
>>> amount = np.asarray([[[5000.0, 0.0]]])
>>> route = np.asarray([[[[[False, True]]], [[[False, True]]], [[[False, True]]]], [[[[True, False]]], [[[False, False]]], [[[False, False]]]], [[[[False, False]]], [[[True, False]]], [[[False, False]]]], [[[[False, False]]], [[[False, False]]], [[[True, False]]]]], dtype=bool)
>>> amount.shape
(1, 1, 2)
>>> Added dimension `x`
>>> # y,x,a,b,c
>>> route.shape
(4, 3, 1, 1, 2)
>>> # Attempt 1: `5000` can flow into y=1, 2 or 3. I expect
>>> # `flows1.sum() == amount.sum()` as it would be without `x`.
>>> # Correct solution would be `[0, 5000, 0, 0]` because material is routed
>>> # to y=1, and is not available for y=2 and y=3 as they are lower
>>> # priority (higher index)
>>> flows1 = np.einsum('abc,yxabc->y', amount, route)
>>> flows1
array([ 0., 5000., 5000., 5000.])
>>> # Attempt 2: try to collapse `x` => not much different, duplication
>>> np.einsum('abc,yabc->y', amount, route.any(1))
array([ 0., 5000., 5000., 5000.])
>>> # This is the flow by `y` and `x`. I'd only expect a `5000` in the
>>> # 2nd row (`[5000., 0., 0.]`) not the others.
>>> np.einsum('abc,yxabc->yx', amount, route)
array([[ 0., 0., 0.],
[5000., 0., 0.],
[ 0., 5000., 0.],
[ 0., 0., 5000.]])
Существует ли какая-либо допустимая операция, которую я могу применить к route
(.all(1)
тоже не работает), чтобы игнорировать размерность x?
Другой пример:
>>> amount2 = np.asarray([[[5000.0, 1000.0]]])
>>> np.einsum('abc,yabc->y', amount2, route.any(1))
array([1000., 5000., 5000., 5000.])
можно интерпретировать как 1000.0
, маршрутизируемый на y=0
(и ни один из других y-пунктов назначения) и 5000.0
совместимый с пунктом назначения y=1
, y=2
и y=3
, но в идеале я бы хотел показать только 5000.0
в y=1
(поскольку это самый низкий индекс и самый высокий приоритет назначения).
Попытка решения
Ниже работает, но не очень глупо.Было бы здорово, если бы цикл можно было удалить.
# Initialise destination
result = np.zeros((route.shape[0]))
# Calculate flow by maintaining all dimensions (this will cause
# double ups because `x` is not part of `amount2`
temp = np.einsum('abc,yxabc->yxabc', amount2, route)
temp_ixs = np.asarray(np.where(temp))
# For each original amount, find the destination (`y`)
for a, b, c in zip(*np.where(amount2)):
# Find where dimensions `abc` are equal in the destination.
# Take the first vector which contains `yxabc` (we get `yx` as result)
ix = np.where((temp_ixs[2:].T == [a, b, c]).all(axis=1))[0][0]
y_ix = temp_ixs.T[ix][0]
# ignored
x_ix = temp_ixs.T[ix][1]
v = amount2[a, b, c]
# build resulting destination
result[y_ix] += v
# result == array([1000., 5000., 0., 0.])
Другими словами, для каждого значения в amount2
я ищу самые низкие индексы yx
в temp
, чтобы значениеможно записать в result[y] = value
(x игнорируется).
>>> temp = np.einsum('abc,yxabc->yx', amount2, route)
>>> temp
# +--- value=1000 at y=0 => result[0] += 1000
# /
array([[1000., 1000., 1000.],
# +--- value=5000 at y=1 => result[1] += 5000
# /
[5000., 0., 0.],
[ 0., 5000., 0.],
[ 0., 0., 5000.]])
>>> result
array([1000., 5000., 0., 0.])
>>> amount2
array([[[5000., 1000.]]])
Еще одна попытка уменьшить размерность route
:
>>> r = route.any(1)
>>> for x in xrange(1, route.shape[0]):
r[x] = r[x] & (r[:x] == False).all(axis=0)
>>> np.einsum('abc,yabc->y', amount2, r)
array([1000., 5000., 0., 0.])
Это по существу сохраняет вышеупомянутый приоритетдано первым измерением route
.Любой массив с более низким приоритетом (с более высоким индексом) не может содержать значение True, если массив с более высоким приоритетом уже имеет значение True в этом субиндексе.Хотя это намного лучше, чем мой явный подход, было бы замечательно, если бы цикл for x in xrange...
мог быть выражен как операция с вектором-пустышкой.