4 horizontale Doppelpräzisionssummen in einem Durchgang mit AVX

Question

Das Problem kann wie folgt beschrieben werden.

Entrada

__m256d a, b, c, d

Ausgabe

__m256d s = {a[0]+a[1]+a[2]+a[3], b[0]+b[1]+b[2]+b[3], 
             c[0]+c[1]+c[2]+c[3], d[0]+d[1]+d[2]+d[3]}

Bisher geleistete Arbeit

Es schien einfach genug: zwei VHADD mit etwas Mischen dazwischen, aber tatsächlich kann die Kombination aller Permutationen, die AVX bietet, nicht genau die Permutation erzeugen, die nötig ist, um dieses Ziel zu erreichen. Lassen Sie mich das erklären:

VHADD x, a, b => x = {a[0]+a[1], b[0]+b[1], a[2]+a[3], b[2]+b[3]}
VHADD y, c, d => y = {c[0]+c[1], d[0]+d[1], c[2]+c[3], d[2]+d[3]}

Könnte ich x und y auf die gleiche Weise permutieren, um Folgendes zu erhalten

x1 = {a[0]+a[1], a[2]+a[3], c[0]+c[1], c[2]+c[3]}
y1 = {b[0]+b[1], b[2]+b[3], d[0]+d[1], d[2]+d[3]}

dann

VHADD s, x1, y1 => s1 = {a[0]+a[1]+a[2]+a[3], b[0]+b[1]+b[2]+b[3], 
                         c[0]+c[1]+c[2]+c[3], d[0]+d[1]+d[2]+d[3]}

Das ist das gewünschte Ergebnis.

Ich muss also nur herausfinden, wie ich die

x,y => {x[0], x[2], y[0], y[2]}, {x[1], x[3], y[1], y[3]}

Leider bin ich zu dem Schluss gekommen, dass dies mit jeder Kombination von VSHUFPD, VBLENDPD, VPERMILPD, VPERM2F128, VUNPCKHPD, VUNPCKLPD nachweislich unmöglich ist. Der springende Punkt ist, dass es unmöglich ist, u[1] und u[2] in einer Instanz u von __m256d zu vertauschen.

Frage

Ist das wirklich eine Sackgasse? Oder habe ich eine Permutationsanweisung übersehen?

Answer 1

VHADD Die Anweisungen sind dazu bestimmt, von regelmäßigen VADD . Der folgende Code sollte Ihnen das Gewünschte liefern:

// {a[0]+a[1], b[0]+b[1], a[2]+a[3], b[2]+b[3]}
__m256d sumab = _mm256_hadd_pd(a, b);
// {c[0]+c[1], d[0]+d[1], c[2]+c[3], d[2]+d[3]}
__m256d sumcd = _mm256_hadd_pd(c, d);

// {a[0]+a[1], b[0]+b[1], c[2]+c[3], d[2]+d[3]}
__m256d blend = _mm256_blend_pd(sumab, sumcd, 0b1100);
// {a[2]+a[3], b[2]+b[3], c[0]+c[1], d[0]+d[1]}
__m256d perm = _mm256_permute2f128_pd(sumab, sumcd, 0x21);

__m256d sum =  _mm256_add_pd(perm, blend);

Dies ergibt das Ergebnis in 5 Anweisungen. Ich hoffe, ich habe die Konstanten richtig verstanden.

Die von Ihnen vorgeschlagene Permutation ist sicherlich möglich, aber sie erfordert mehrere Anweisungen. Es tut mir leid, dass ich diesen Teil Ihrer Frage nicht beantworte.

Bearbeiten: Ich konnte nicht widerstehen, hier ist die vollständige Permutation. (Auch hier habe ich mein Bestes getan, um die Konstanten richtig hinzubekommen.) Sie können sehen, dass das Vertauschen u[1] y u[2] ist möglich, erfordert nur ein wenig Arbeit. Das Überschreiten der 128-Bit-Grenze ist bei der ersten Generation schwierig. AVX. Ich möchte auch sagen, dass VADD ist vorzuziehen gegenüber VHADD porque VADD hat den doppelten Durchsatz, obwohl es die gleiche Anzahl von Additionen durchführt.

// {x[0],x[1],x[2],x[3]}
__m256d x;

// {x[1],x[0],x[3],x[2]}
__m256d xswap = _mm256_permute_pd(x, 0b0101);

// {x[3],x[2],x[1],x[0]}
__m256d xflip128 = _mm256_permute2f128_pd(xswap, xswap, 0x01);

// {x[0],x[2],x[1],x[3]} -- not imposssible to swap x[1] and x[2]
__m256d xblend = _mm256_blend_pd(x, xflip128, 0b0110);

// repeat the same for y
// {y[0],y[2],y[1],y[3]}
__m256d yblend;

// {x[0],x[2],y[0],y[2]}
__m256d x02y02 = _mm256_permute2f128_pd(xblend, yblend, 0x20);

// {x[1],x[3],y[1],y[3]}
__m256d x13y13 = _mm256_permute2f128_pd(xblend, yblend, 0x31);

Answer 2

Mir ist keine Anweisung bekannt, die diese Art von Permutation ermöglicht. AVX-Befehle arbeiten in der Regel so, dass die oberen und unteren 128 Bits des Registers einigermaßen unabhängig sind; es gibt nicht viele Möglichkeiten, Werte aus den beiden Hälften zu vermischen. Die beste Implementierung, die ich mir vorstellen kann, würde auf der Antwort auf diese Frage :

__m128d horizontal_add_pd(__m256d x1, __m256d x2)
{
    // calculate 4 two-element horizontal sums:
    // lower 64 bits contain x1[0] + x1[1]
    // next 64 bits contain x2[0] + x1[1]
    // next 64 bits contain x1[2] + x1[3]
    // next 64 bits contain x2[2] + x2[3]
    __m256d sum = _mm256_hadd_pd(x1, x2);
    // extract upper 128 bits of result
    __m128d sum_high = _mm256_extractf128_pd(sum1, 1);
    // add upper 128 bits of sum to its lower 128 bits
    __m128d result = _mm_add_pd(sum_high, (__m128d) sum);
    // lower 64 bits of result contain the sum of x1[0], x1[1], x1[2], x1[3]
    // upper 64 bits of result contain the sum of x2[0], x2[1], x2[2], x2[3]
    return result;
}

__m256d a, b, c, d;
__m128d res1 = horizontal_add_pd(a, b);
__m128d res2 = horizontal_add_pd(c, d);
// At this point:
//     res1 contains a's horizontal sum in bits 0-63
//     res1 contains b's horizontal sum in bits 64-127
//     res2 contains c's horizontal sum in bits 0-63
//     res2 contains d's horizontal sum in bits 64-127
// cast res1 to a __m256d, then insert res2 into the upper 128 bits of the result
__m256d sum = _mm256_insertf128_pd(_mm256_castpd128_pd256(res1), res2, 1);
// At this point:
//     sum contains a's horizontal sum in bits 0-63
//     sum contains b's horizontal sum in bits 64-127
//     sum contains c's horizontal sum in bits 128-191
//     sum contains d's horizontal sum in bits 192-255

Das sollte das sein, was Sie wollen. Das obige sollte in insgesamt 7 Anweisungen machbar sein (der Cast sollte nicht wirklich etwas tun; es ist nur ein Hinweis an den Compiler, die Art und Weise zu ändern, wie er den Wert in res1 ), unter der Annahme, dass die kurze horizontal_add_pd() Funktion von Ihrem Compiler eingebettet werden kann und Sie genügend Register zur Verfügung haben.