Почему?Дорого ли _mm_setr_pd?
В некотором смысле;это займет хотя бы случайное перемешивание.Что еще более важно в этом случае, вычисление каждого скалярного операнда стоит дорого, и, как показывает ответ @spectras, gcc по крайней мере не может автоматически векторизовать это в paddd
/ cvtdq2pd
.Вместо этого он заново вычисляет каждый операнд из скалярного целого числа, выполняя преобразование int
-> double
отдельно, а затем перемешивает их вместе.
Это функция, которую я пытаюсь оптимизировать:
Вы просто заполняете массив линейной функцией .Вы повторяете каждый раз внутри цикла.Это позволяет избежать зависимости, переносимой в цикле, от чего угодно, кроме целочисленного счетчика цикла, но вы сталкиваетесь с узкими местами пропускной способности из-за такой большой работы внутри цикла.
т.е. вы вычисляете a[i] = c + i*scale
отдельно для каждого шага.Но вместо вы можете уменьшить его до a[i+n] = a[i] + (n*scale)
.Таким образом, у вас есть только одна addpd
инструкция на вектор результатов.
Это приведет к некоторой ошибке округления, которая накапливается по сравнению с повторным вычислением с нуля, но double
, вероятно, излишне для вас.в любом случае.
Это также происходит за счет введения последовательной зависимости от добавления FP вместо целого числа.Но в вашей «оптимизированной» версии уже есть цепочка зависимостей добавления FP, переносимая в цикле, которая использует внутри цикла sampleIndex_acc = _mm_add_pd(sampleIndex_acc, sampleIndex_add);
, используя FP + = 2.0 вместо повторного преобразования из целого числа.
Итаквам нужно будет развернуть несколько векторов, чтобы скрыть задержку FP, и сохранить как минимум 3 или 4 дополнения FP в полете одновременно .(Haswell: задержка 3 цикла, один на тактовую пропускную способность. Skylake: задержка 4 цикла, 2 на тактовую пропускную способность.) См. Также Почему mulss занимает только 3 цикла на Haswell, в отличие от таблиц инструкций Агнера? для получения дополнительной информациио развертывании с несколькими аккумуляторами для аналогичной проблемы с зависимостями, переносимыми циклом (точечный продукт).
void Process(int voiceIndex, int blockSize) {
double *pC = c[voiceIndex];
double val0 = start + step * delta;
double deltaValue = rate * delta;
__m128d vdelta2 = _mm_set1_pd(2 * deltaValue);
__m128d vdelta4 = _mm_add_pd(vdelta2, vdelta2);
__m128d v0 = _mm_setr_pd(val0, val0 + deltaValue);
__m128d v1 = _mm_add_pd(v0, vdelta2);
__m128d v2 = _mm_add_pd(v0, vdelta4);
__m128d v3 = _mm_add_pd(v1, vdelta4);
__m128d vdelta8 = _mm_mul_pd(vdelta2, _mm_set1_pd(4.0));
double *endp = pC + blocksize - 7; // stop if there's only room for 7 or fewer doubles
// or use -8 and have your cleanup handle lengths of 1..8
// since the inner loop always calculates results for next iteration
for (; pC < endp ; pC += 8) {
_mm_store_pd(pC, v0);
v0 = _mm_add_pd(v0, vdelta8);
_mm_store_pd(pC+2, v1);
v1 = _mm_add_pd(v1, vdelta8);
_mm_store_pd(pC+4, v2);
v2 = _mm_add_pd(v2, vdelta8);
_mm_store_pd(pC+6, v3);
v3 = _mm_add_pd(v3, vdelta8);
}
// if (blocksize % 8 != 0) ... store final vectors
}
Выбор, добавлять или умножать при создании vdelta4
/ vdelta8
, не оченьзначительное;Я пытался избежать слишком длинной цепочки зависимостей, прежде чем появятся первые магазины.Поскольку от 1041 * до v3
также необходимо вычислить, казалось, что имеет смысл создать vdelta4
вместо того, чтобы просто создавать цепочку v2 = v1+vdelta2
.Возможно, было бы лучше создать vdelta4
с умножением от 4.0*delta
и удвоить его, чтобы получить vdelta8
.Это может иметь отношение к очень маленькому размеру блока, особенно если вы кешируете свой код, генерируя только небольшие куски этого массива по мере необходимости, прямо перед его чтением.
В любом случае, это компилируется в очень эффективныйвнутренний цикл с gcc и MSVC ( в проводнике компилятора Godbolt ).
;; MSVC -O2
$LL4@Process: ; do {
movups XMMWORD PTR [rax], xmm5
movups XMMWORD PTR [rax+16], xmm0
movups XMMWORD PTR [rax+32], xmm1
movups XMMWORD PTR [rax+48], xmm2
add rax, 64 ; 00000040H
addpd xmm5, xmm3 ; v0 += vdelta8
addpd xmm0, xmm3 ; v1 += vdelta8
addpd xmm1, xmm3 ; v2 += vdelta8
addpd xmm2, xmm3 ; v3 += vdelta8
cmp rax, rcx
jb SHORT $LL4@Process ; }while(pC < endp)
Это имеет 4 отдельных цепочки зависимостей, через xmm0, 1, 2 и 5. Так что достаточно инструкции-уровень параллелизма для сохранения 4 addpd
инструкций в полете.Этого более чем достаточно для Haswell, но половина того, что может выдержать Skylake.
Тем не менее, с пропускной способностью магазина 1 вектор на такт, более 1 addpd
на такт бесполезна. Теоретически это может быть около 16 байтов за такт, и насыщать пропускную способность хранилища. т.е. 1 вектор / 2 double
с на такт.
AVX с более широкими векторами (4 double
s) все еще может идти с 1 вектором за такт в Haswell и позже, то есть 32 байта за такт.(Предполагая, что выходной массив горячий в кеше L1d или, возможно, даже в L2.)
Еще лучше: вообще не хранить эти данные в памяти;сгенерируйте на лету.
Создайте его на лету, когда это необходимо, если код, потребляющий его, только читает его несколько раз, а также векторизует вручную.