Неинтуитивно понятный параметр сравнения между `matrix * vector`,` matrix '* vector` и `copy (matrix') * vector`

Question

Вопрос от Джулия Дискурс

Я использую Джулию 1.2.Это мой тест:

a = rand(1000, 1000)
b = a'
c = copy(b)

@btime a * x setup=(x=rand(1000)) # 114.757 μs
@btime b * x setup=(x=rand(1000)) # 94.179 μs
@btime c * x setup=(x=rand(1000)) # 110.325 μs

Я ожидал, что a и c будут как минимум не медленнее.

После проверки stdlib / LinearAlgebra / src / matmul.jl , оказывается, что Юлия передает b.parent (то есть a ) в BLAS.gemv ,не b , а вместо этого переключает dgemv_ LAPACK в другой и, по-видимому, более быстрый режим.

Правильно ли я считаю, что ускорение происходит из-за того, что памятьвыровнен более благоприятным образом для всего, что dgemv_ делает, когда он находится в режиме trans = T ?Если так, то я предполагаю, что это невозможно, кроме возможного упоминания о гоче в документах.Хотя, если мое предположение неверно, нужно ли что-то делать с этим?

Answer 1

Ответ от @stevengj в той же ветке Discourse:

Правильно ли я считаю, что ускорение происходит из-за того, что память выровнена более благоприятным образом для всего, что угодно dgemv_ делает, когда он находится в режиме trans = T ?

Закрыть.Это связано с памятью, но речь идет о locality , а не о выравнивании.Основная вещь, которую нужно понять, это то, что более эффективно получать доступ к последовательным (или, по крайней мере, поблизости) данным из памяти, чем к данным, которые разделены, из-за существования строк кэша .(Последовательный доступ также имеет некоторые преимущества при использовании SIMD инструкций.)

Юлия хранит матрицы в мажорном порядке столбцов , так что столбцы непрерывны в памяти.Следовательно, когда вы умножаете транспонированную матрицу (которая не была скопирована) на вектор, она может вычислить ее как произведение точек смежного столбца (= транспонированная строка) на смежный вектор, который имеет хороший пространственная локальность и, следовательно, эффективно использует строки кэша.

Для умножения не -транспонированной матрицы на вектор, напротив, вы берете точечные произведения несмежные строки матрицы с вектором, и сложнее эффективно использовать строки кэша.Я полагаю, что для улучшения пространственной локальности оптимизированный BLAS, такой как OpenBLAS, фактически вычисляет точечные произведения нескольких строк за раз («блок») с вектором - поэтому он всего на 10% медленнее и не намного хуже.(Фактически, даже транспонированный случай может сделать некоторую блокировку, чтобы сохранить вектор в кэше.)