Создание матрицы из функции с несколькими массивами без циклов for - PullRequest
1 голос
/ 03 апреля 2020

У меня есть следующие пять массивов:

t = np.linspace(0,100,100)
Q = t**2

x = np.linspace(10,20,20)
y = np.linspace(5,8,5)
z = np.linspace(100,125,10)

Я также определяю функцию следующим образом:

def f(t, x, y, z):
    return (1+x)**3.2 * y**x * np.cos(y**x)*(t-z)

Теперь для каждого i-го элемента x, j- й элемент y и k-й элемент z, я хочу определить значение, которое f(t, x, y, z) принимает при всех значениях t. Я также хочу взять «внутренний продукт» с каждым из этих результатов с массивом Q. Итак, в итоге я хочу получить трехмерную матрицу, которая хранит массивы размером len(t), где каждый элемент (i, j, k) этой матрицы равен np.sum( Q * f(t, x[i], y[j], z[k])). Все это может быть достигнуто с помощью a для l oop следующим образом:

result = np.zeros(len(x)*len(y)*len(z)).reshape((len(x),len(y), len(z)))
for i in range(len(x)):
    for j in range(len(y)):
        for k in range(len(z)):
            intermediate_result = f(t, x[i], y[j], z[k]) #This is an array of length len(t)
            result[i][j][k] = (np.nansum(intermediate_result * Q))**2

Проблема состоит в том, что эти циклы for требуют очень много времени для вычисления больших массивов x, y и z за пределами этого упрощенного примера. поэтому я ищу способ уменьшить время вычислений. Есть ли эффективный способ сделать это?

1 Ответ

2 голосов
/ 04 апреля 2020

Вам не нужно все oop. вы можете позволить numpy осуществлять трансляцию, изменяя параметры переменных:

result = f(t,x[:,None,None,None],y[:,None,None],z[:,None])

result = np.sum(Q*result,axis=3)

обратите внимание, что вы можете сделать это внутри функции, используя np.ix_

def f(t, x, y, z):
    x,y,z,t = np.ix_(x,y,z,t)
    return (1+x)**3.2 * y**x * np.cos(y**x)*(t-z)

Функция вернет 4D матрица, которую вы можете затем обработать с помощью np.sum:

result = np.sum(Q*f(t,x,y,z),axis=3)

Это будет работать заметно быстрее, чем циклы (более чем в 10 раз быстрее для этого примера 100 x 20x5x10)

[РЕДАКТИРОВАТЬ] Учитывая, что 4D матрица занимает много памяти, вы в конечном итоге попадете в другое узкое место, когда размеры станут больше. Вы можете обойти это, ограничив матрицу тремя измерениями и вручную добавив измерение t в al oop. Это гарантирует, что вы никогда не будете использовать размер получаемой трехмерной матрицы более чем в два раза:

def f2(t, x, y, z):
    x,y,z = np.ix_(x,y,z)
    return (1+x)**3.2 * y**x * np.cos(y**x)*(t-z)

result = np.zeros((x.size,y.size,z.size))
for tn,Qn in zip(t,Q):
    result += f2(tn,x,y,z)*Qn

При таком подходе мне удалось получить результат за 19 секунд при t.size = 300 и x, y, z при 256х256х256. Матричный подход 4D разбил оболочку IDLE (возможно, из-за переполнения памяти), и у меня не хватило терпения дождаться окончания цикла оригинала до 1027 *.

. Обратите внимание, что этот последний подход может привести к незначительным различиям около 16-го числа git мантиссы по сравнению с исходным решением из-за порядка, в котором np.sum () добавляет значения

...