Понимание std :: transform и как с этим справиться

Я пытаюсь понять ориентированный на данные дизайн простой, специфической проблемы. Заранее извиняюсь перед людьми, ориентированными на данные, если я делаю что-то очень глупое, но мне трудно понять, почему и где мои рассуждения терпят неудачу.

Предположим, что у меня есть простая операция, то есть, float_t result = int_t(lhs) / int_t(rhs), Если я храню все переменные в соответствующих контейнерах, например, std::vector<float_t> а также std::vector<int_t>и я использую std::transformЯ получаю правильный результат. Затем для конкретного примера, где using float_t = float а также using int_t = int16_tЯ предполагаю, что упаковка этих переменных внутри structна 64-битной архитектуре и собирая их внутри контейнера должен дать лучшую производительность.

Я считаю, что struct составляет 64-битный объект и один доступ к памяти struct даст мне все переменные, которые мне нужны. С другой стороны, когда все эти переменные собраны в разных контейнерах, мне потребуется три разных доступа к памяти, чтобы получить необходимую информацию. Ниже описано, как я настраиваю среду:

#include <algorithm>
#include <chrono>
#include <cstdint>
#include <iostream>
#include <vector>

using namespace std::chrono;

template <class float_t, class int_t> struct Packed {
float_t sinvl;
int_t s, l;
Packed() = default;
Packed(float_t sinvl, int_t s, int_t l) : sinvl{sinvl}, s{s}, l{l} {}
void comp() { sinvl = float_t(l) / s; }
};

using my_float = float;
using my_int = int16_t;

int main(int argc, char *argv[]) {
constexpr uint32_t M{100};
for (auto N : {1000, 10000, 100000}) {
double t1{0}, t2{0};
for (uint32_t m = 0; m < M; m++) {
std::vector<my_float> sinvl(N, 0.0);
std::vector<my_int> s(N, 3), l(N, 2);
std::vector<Packed<my_float, my_int>> p1(
N, Packed<my_float, my_int>(0.0, 3, 2));

// benchmark unpacked
auto tstart = high_resolution_clock::now();
std::transform(l.cbegin(), l.cend(), s.cbegin(), sinvl.begin(),
std::divides<my_float>{}); // 3 different memory accesses
auto tend = high_resolution_clock::now();
t1 += duration_cast<microseconds>(tend - tstart).count();

if (m == M - 1)
std::cout << "sinvl[0]: " << sinvl[0] << '\n';

// benchmark packed
tstart = high_resolution_clock::now();
for (auto &elem : p1) // 1 memory access
elem.comp();
tend = high_resolution_clock::now();
t2 += duration_cast<microseconds>(tend - tstart).count();

if (m == M - 1)
std::cout << "p1[0].sinvl: " << p1[0].sinvl << '\n';
}
std::cout << "N = " << N << ", unpacked: " << (t1 / M) << " us.\n";
std::cout << "N = " << N << ", packed: " << (t2 / M) << " us.\n";
}
return 0;
}

Скомпилированный код с g++ -O3 дает, на моей машине,

sinvl[0]: 0.666667
p1[0].sinvl: 0.666667
N = 1000, unpacked: 0 us.
N = 1000, packed: 1 us.
sinvl[0]: 0.666667
p1[0].sinvl: 0.666667
N = 10000, unpacked: 5.06 us.
N = 10000, packed: 12.97 us.
sinvl[0]: 0.666667
p1[0].sinvl: 0.666667
N = 100000, unpacked: 52.31 us.
N = 100000, packed: 124.49 us.

В принципе, std::transform бьет упакованный доступ 2.5x, Я был бы признателен, если бы вы помогли мне понять поведение. Является ли результат из-за

  1. я не правильно понимаю принципы ориентированного на данные дизайна, или,
  2. какой-нибудь артефакт этого очень простого примера, такой как области памяти, расположенные очень близко друг к другу и каким-то образом оптимизируемые компилятором очень эффективно?

Наконец, есть ли способ победить std::transform в этом примере, или, это просто достаточно хорошо, чтобы быть готовым решением? Я не являюсь экспертом ни в оптимизации компиляторов, ни в ориентированном на данные проектировании, и, таким образом, я не мог ответить на этот вопрос сам.

Спасибо!

РЕДАКТИРОВАТЬ. Я изменил способ тестирования обоих методов согласно предложению @ bolov в комментариях.

Теперь код выглядит так:

#include <algorithm>
#include <chrono>
#include <cstdint>
#include <iostream>
#include <vector>

using namespace std::chrono;

template <class float_t, class int_t> struct Packed {
float_t sinvl;
int_t s, l;
Packed() = default;
Packed(float_t sinvl, int_t s, int_t l) : sinvl{sinvl}, s{s}, l{l} {}
void comp() { sinvl = float_t(l) / s; }
};

using my_float = float;
using my_int = int16_t;

int main(int argc, char *argv[]) {
uint32_t N{1000};
double t{0};

if (argc == 2)
N = std::stoul(argv[1]);

#ifndef _M_PACKED
std::vector<my_float> sinvl(N, 0.0);
std::vector<my_int> s(N, 3), l(N, 2);

// benchmark unpacked
auto tstart = high_resolution_clock::now();
std::transform(l.cbegin(), l.cend(), s.cbegin(), sinvl.begin(),
std::divides<my_float>{}); // 3 different memory accesses
auto tend = high_resolution_clock::now();
t += duration_cast<microseconds>(tend - tstart).count();

std::cout << "sinvl[0]: " << sinvl[0] << '\n';
std::cout << "N = " << N << ", unpacked: " << t << " us.\n";
#else
std::vector<Packed<my_float, my_int>> p1(N,
Packed<my_float, my_int>(0.0, 3, 2));
// benchmark packed
auto tstart = high_resolution_clock::now();
for (auto &elem : p1) // 1 memory access
elem.comp();
auto tend = high_resolution_clock::now();
t += duration_cast<microseconds>(tend - tstart).count();

std::cout << "p1[0].sinvl: " << p1[0].sinvl << '\n';
std::cout << "N = " << N << ", packed: " << t << " us.\n";
#endif

return 0;
}

с соответствующим сценарием оболочки (рыбы)

g++ -Wall -std=c++11 -O3 transform.cpp -o transform_unpacked.out
g++ -Wall -std=c++11 -O3 transform.cpp -o transform_packed.out -D_M_PACKED
for N in 1000 10000 100000
echo "Testing unpacked for N = $N"./transform_unpacked.out $N
./transform_unpacked.out $N
./transform_unpacked.out $N
echo "Testing packed for N = $N"./transform_packed.out $N
./transform_packed.out $N
./transform_packed.out $N
end

что дает следующее:

Testing unpacked for N = 1000
sinvl[0]: 0.666667
N = 1000, unpacked: 0 us.
sinvl[0]: 0.666667
N = 1000, unpacked: 0 us.
sinvl[0]: 0.666667
N = 1000, unpacked: 0 us.
Testing packed for N = 1000
p1[0].sinvl: 0.666667
N = 1000, packed: 1 us.
p1[0].sinvl: 0.666667
N = 1000, packed: 1 us.
p1[0].sinvl: 0.666667
N = 1000, packed: 1 us.
Testing unpacked for N = 10000
sinvl[0]: 0.666667
N = 10000, unpacked: 5 us.
sinvl[0]: 0.666667
N = 10000, unpacked: 5 us.
sinvl[0]: 0.666667
N = 10000, unpacked: 5 us.
Testing packed for N = 10000
p1[0].sinvl: 0.666667
N = 10000, packed: 17 us.
p1[0].sinvl: 0.666667
N = 10000, packed: 13 us.
p1[0].sinvl: 0.666667
N = 10000, packed: 13 us.
Testing unpacked for N = 100000
sinvl[0]: 0.666667
N = 100000, unpacked: 64 us.
sinvl[0]: 0.666667
N = 100000, unpacked: 66 us.
sinvl[0]: 0.666667
N = 100000, unpacked: 66 us.
Testing packed for N = 100000
p1[0].sinvl: 0.666667
N = 100000, packed: 180 us.
p1[0].sinvl: 0.666667
N = 100000, packed: 198 us.
p1[0].sinvl: 0.666667
N = 100000, packed: 177 us.

Надеюсь, я правильно понял метод тестирования. Тем не менее, разница составляет 2-3 раза.

10

Решение

Вот скомпилированный цикл std::transform дело:

  400fd0:       f3 41 0f 7e 04 47       movq   xmm0,QWORD PTR [r15+rax*2]
400fd6:       66 0f 61 c0             punpcklwd xmm0,xmm0
400fda:       66 0f 72 e0 10          psrad  xmm0,0x10
400fdf:       0f 5b c0                cvtdq2ps xmm0,xmm0
400fe2:       f3 0f 7e 0c 43          movq   xmm1,QWORD PTR [rbx+rax*2]
400fe7:       66 0f 61 c9             punpcklwd xmm1,xmm1
400feb:       66 0f 72 e1 10          psrad  xmm1,0x10
400ff0:       0f 5b c9                cvtdq2ps xmm1,xmm1
400ff3:       0f 5e c1                divps  xmm0,xmm1
400ff6:       41 0f 11 04 80          movups XMMWORD PTR [r8+rax*4],xmm0
400ffb:       f3 41 0f 7e 44 47 08    movq   xmm0,QWORD PTR [r15+rax*2+0x8]
401002:       66 0f 61 c0             punpcklwd xmm0,xmm0
401006:       66 0f 72 e0 10          psrad  xmm0,0x10
40100b:       0f 5b c0                cvtdq2ps xmm0,xmm0
40100e:       f3 0f 7e 4c 43 08       movq   xmm1,QWORD PTR [rbx+rax*2+0x8]
401014:       66 0f 61 c9             punpcklwd xmm1,xmm1
401018:       66 0f 72 e1 10          psrad  xmm1,0x10
40101d:       0f 5b c9                cvtdq2ps xmm1,xmm1
401020:       0f 5e c1                divps  xmm0,xmm1
401023:       41 0f 11 44 80 10       movups XMMWORD PTR [r8+rax*4+0x10],xmm0
401029:       48 83 c0 08             add    rax,0x8
40102d:       48 83 c1 02             add    rcx,0x2
401031:       75 9d                   jne    400fd0 <main+0x570>

В каждом цикле цикла обрабатывается 8 элементов (есть два divps инструкции, каждый делает 4 деления).

Вот другой случай:

  401190:       f3 0f 6f 42 04          movdqu xmm0,XMMWORD PTR [rdx+0x4]
401195:       f3 0f 6f 4a 14          movdqu xmm1,XMMWORD PTR [rdx+0x14]
40119a:       66 0f 70 c9 e8          pshufd xmm1,xmm1,0xe8
40119f:       66 0f 70 c0 e8          pshufd xmm0,xmm0,0xe8
4011a4:       f2 0f 70 d0 e8          pshuflw xmm2,xmm0,0xe8
4011a9:       66 0f 6c c1             punpcklqdq xmm0,xmm1
4011ad:       66 0f 72 e0 10          psrad  xmm0,0x10
4011b2:       0f 5b c0                cvtdq2ps xmm0,xmm0
4011b5:       f2 0f 70 c9 e8          pshuflw xmm1,xmm1,0xe8
4011ba:       66 0f 62 d1             punpckldq xmm2,xmm1
4011be:       66 0f 61 ca             punpcklwd xmm1,xmm2
4011c2:       66 0f 72 e1 10          psrad  xmm1,0x10
4011c7:       0f 5b c9                cvtdq2ps xmm1,xmm1
4011ca:       0f 5e c1                divps  xmm0,xmm1
4011cd:       f3 0f 11 02             movss  DWORD PTR [rdx],xmm0
4011d1:       0f 28 c8                movaps xmm1,xmm0
4011d4:       0f c6 c9 e5             shufps xmm1,xmm1,0xe5
4011d8:       f3 0f 11 4a 08          movss  DWORD PTR [rdx+0x8],xmm1
4011dd:       0f 28 c8                movaps xmm1,xmm0
4011e0:       0f 12 c9                movhlps xmm1,xmm1
4011e3:       f3 0f 11 4a 10          movss  DWORD PTR [rdx+0x10],xmm1
4011e8:       0f c6 c0 e7             shufps xmm0,xmm0,0xe7
4011ec:       f3 0f 11 42 18          movss  DWORD PTR [rdx+0x18],xmm0
4011f1:       48 83 c2 20             add    rdx,0x20
4011f5:       48 83 c1 fc             add    rcx,0xfffffffffffffffc
4011f9:       75 95                   jne    401190 <main+0x730>

В каждом цикле цикла обрабатывается 4 элемента (есть один divps инструкция).

В первом случае данные находятся в хорошем формате, SIMD-инструкции могут работать с ними (почти) без какого-либо перемещения данных, и результат может быть легко записан (4 результата записываются с помощью одной инструкции).

Во втором случае, однако, это не так. Компилятор должен был выполнить много операций перемещения данных (случайного перемещения), и каждый результат записывается с отдельной инструкцией. Поэтому ввод / вывод не в формате SIMD.

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

4

Другие решения

[…] и собирая их в контейнере, вы получите лучший результат
спектакль.

Я думаю, что ваша интуиция немного отсталая для последовательного доступа. Если учесть последовательные циклы на входах нетривиального размера, SoA почти всегда будет работать быстрее или быстрее, чем AoS для последовательный доступ.

Во многих случаях SoA на самом деле приведет к меньшему количеству пропусков кэша, чем AoS, потому что он избегает заполнения структуры (не отвечает требованиям выравнивания чередования неоднородных полей), и вы можете избежать загрузки холодных полей для данного цикла (например, : симулятор физики может избежать загрузки color поле частицы, используемое только для рендеринга, при применении физики). Ваш случай не выигрывает ни от одного из них, но об этом следует помнить.

Где AoS стремится к совершенству для произвольный доступ. В этом случае, если вы загрузите, скажем, 4096-й элемент, вам может потребоваться только одна строка кэша с AoS, чтобы получить все три поля. Если вы используете SoA, вам потребуется 3, и он может загрузить много данных для соседних элементов (данных для элемента 4097, 4098, 4099 и т. Д.), Ни один из которых не используется до выселения из-за характера произвольного доступа шаблон доступа к памяти. Но для последовательного доступа все такие соседние данные, как правило, будут использоваться даже с представителем SoA до выселения.

Таким образом, с точки зрения промахов кэша в тех случаях, когда вы последовательно зацикливаетесь на нетривиальных размерах входных данных, SoA обычно либо связывает, либо выигрывает.

Но, кроме того, в таких последовательных схемах доступа, где все соседние данные будут использоваться до выселения, при рассмотрении динамики загрузки данных из кэша в регистр SIMD SoA предоставляет однородные данные. Он может загружать память в виде, скажем, float, float, float, float, ... а также int16, int16, int16, int16, ...не float, int16, int16, float, int16, int16 ... с вертикальными / симметричными операциями, выполняемыми через регистры SIMD. Такие случаи, как правило, предлагают гораздо больше возможностей для оптимизации как для человека, так и для компилятора, и это, вероятно, когда ваш конкретный случай приносит наибольшую пользу, как указывает Геза.

Наконец, есть ли способ превзойти std :: transform в этом примере, или
это просто достаточно хорошо, чтобы быть готовым решением?

Я думаю, пытаясь победить std::transform хотя выполнение тех же требований является неправильным подходом, вы можете уменьшить некоторые приросты производительности, изменив их немного. Например, вместо std::dividesВы можете подготовить данные заранее, чтобы превратить их в умножение. Самая большая вещь, которую я бы искал в вашем случае, это посмотреть, может ли код выполняться параллельно с более быстрым представителем SoA («распакованный»). Возможно, вы сможете обработать заданный диапазон индексов каждого контейнера SoA в каждом отдельном потоке, по-прежнему используя std::transform в каждом.

1