Реализация NumPy тензор медленнее, чем цикл - PullRequest
Новая
       12

Реализация NumPy тензор медленнее, чем цикл

1 голос
/ 21 марта 2019

У меня есть две функции, которые вычисляют одну и ту же метрику. Один из них в конечном итоге использует понимание списка для циклического вычисления, другой использует только операции с тензорными числами. Функции принимают массив (N, 3), где N - количество точек в трехмерном пространстве. Когда N <~ 3000, тензорная функция быстрее, когда N> ~ 3000, понимание списка происходит быстрее. Кажется, что обе имеют линейную временную сложность с точки зрения N, то есть две линии время-N пересекаются при N = ~ 3000. 

def approximate_area_loop(section, num_area_divisions):        
    n_a_d = num_area_divisions
    interp_vectors = get_section_interp_(section)

    a1 = section[:-1]
    b1 = section[1:]
    a2 = interp_vectors[:-1]
    b2 = interp_vectors[1:]

    c = lambda u: (1 - u) * a1 + u * a2
    d = lambda u: (1 - u) * b1 + u * b2
    x = lambda u, v: (1 - v) * c(u) + v * d(u)

    area = np.sum([np.linalg.norm(np.cross((x((i + 1)/n_a_d, j/n_a_d) - x(i/n_a_d, j/n_a_d)),\
                                           (x(i/n_a_d, (j +1)/n_a_d) - x(i/n_a_d, j/n_a_d))), axis = 1)\
                   for i in range(n_a_d) for j in range(n_a_d)])

    Dt = section[-1, 0] - section[0, 0]
    return area, Dt

def approximate_area_tensor(section, num_area_divisions):
    divisors = np.linspace(0, 1, num_area_divisions + 1)
    interp_vectors = get_section_interp_(section)
    a1 = section[:-1]
    b1 = section[1:]
    a2 = interp_vectors[:-1]
    b2 = interp_vectors[1:]
    c = np.multiply.outer(a1, (1 - divisors)) + np.multiply.outer(a2, divisors) # c_areas_vecs_divs
    d = np.multiply.outer(b1, (1 - divisors)) + np.multiply.outer(b2, divisors) # d_areas_vecs_divs
    x = np.multiply.outer(c, (1 - divisors)) + np.multiply.outer(d, divisors) # x_areas_vecs_Divs_divs
    u = x[:, :, 1:, :-1] - x[:, :, :-1, :-1] # u_areas_vecs_Divs_divs
    v = x[:, :, :-1, 1:] - x[:, :, :-1, :-1] # v_areas_vecs_Divs_divs
    sub_area_norm_vecs = np.cross(u, v, axis = 1) # areas_crosses_Divs_divs
    sub_areas = np.linalg.norm(sub_area_norm_vecs, axis = 1) # areas_Divs_divs (values are now sub areas)
    area = np.sum(sub_areas)
    Dt = section[-1, 0] - section[0, 0]
    return area, Dt

Почему версия со списком работает быстрее при больших N? Неужели тензорная версия должна быть быстрее? Мне интересно, связано ли это с размером вычислений, означающим, что он слишком велик, чтобы его можно было сделать в кеше? Пожалуйста, спросите, если я не включил достаточно информации, я бы очень хотел, чтобы до конца это.

1 Ответ

0 голосов
/ 26 марта 2019

Узкое место в полностью векторизованной функции действительно было в np.linalg.norm, как было предложено в комментарии @hpauljs. Норма использовалась только для получения величины всех векторов, содержащихся в оси 1. Гораздо более простой и быстрый метод заключался в следующем: 

sub_areas = np.sqrt((sub_area_norm_vecs*sub_area_norm_vecs).sum(axis = 1))

Это дает точно такие же результаты и ускоряет код в 25 раз быстрее, чем реализация цикла (даже если цикл также не использует linalg.norm).

...