Почему это умножение SIMD не быстрее, чем умножение не SIMD?
Давайте предположим, что у нас есть функция, которая умножает два массива на 1000000 удваивается каждый. В C/C++ функция выглядит так:
void mul_c(double* a, double* b)
{
for (int i = 0; i != 1000000; ++i)
{
a[i] = a[i] * b[i];
}
}
Компилятор производит следующую сборку с -O2:
mul_c(double*, double*):
xor eax, eax
.L2:
movsd xmm0, QWORD PTR [rdi+rax]
mulsd xmm0, QWORD PTR [rsi+rax]
movsd QWORD PTR [rdi+rax], xmm0
add rax, 8
cmp rax, 8000000
jne .L2
rep ret
Из вышеприведенной сборки кажется, что компилятор использует SIMD-инструкции, но он умножает только одну двойную на каждую итерацию. Поэтому я решил написать ту же функцию во встроенной сборке, где я полностью использую xmm0 зарегистрируйте и умножьте два двойных за один раз:
void mul_asm(double* a, double* b)
{
asm volatile
(
".intel_syntax noprefix \n\t"
"xor rax, rax \n\t"
"0: \n\t"
"movupd xmm0, xmmword ptr [rdi+rax] \n\t"
"mulpd xmm0, xmmword ptr [rsi+rax] \n\t"
"movupd xmmword ptr [rdi+rax], xmm0 \n\t"
"add rax, 16 \n\t"
"cmp rax, 8000000 \n\t"
"jne 0b \n\t"
".att_syntax noprefix \n\t"
:
: "D" (a), "S" (b)
: "memory", "cc"
);
}
После измерения времени выполнения по отдельности для обеих этих функций кажется, что для выполнения каждой из них требуется 1 мс:
> gcc -O2 main.cpp
> ./a.out < input
mul_c: 1 ms
mul_asm: 1 ms
[a lot of doubles...]
Я ожидал, что реализация SIMD будет как минимум вдвое быстрее (0 мс), так как есть только половина количества операций умножения / памяти.
Итак, мой вопрос: почему реализация SIMD не быстрее, чем обычная реализация C/C++, когда реализация SIMD выполняет только половину операций умножения / памяти?
Вот полная программа:
#include <stdio.h>
#include <stdlib.h>
#include <sys/time.h>
void mul_c(double* a, double* b)
{
for (int i = 0; i != 1000000; ++i)
{
a[i] = a[i] * b[i];
}
}
void mul_asm(double* a, double* b)
{
asm volatile
(
".intel_syntax noprefix \n\t"
"xor rax, rax \n\t"
"0: \n\t"
"movupd xmm0, xmmword ptr [rdi+rax] \n\t"
"mulpd xmm0, xmmword ptr [rsi+rax] \n\t"
"movupd xmmword ptr [rdi+rax], xmm0 \n\t"
"add rax, 16 \n\t"
"cmp rax, 8000000 \n\t"
"jne 0b \n\t"
".att_syntax noprefix \n\t"
:
: "D" (a), "S" (b)
: "memory", "cc"
);
}
int main()
{
struct timeval t1;
struct timeval t2;
unsigned long long time;
double* a = (double*)malloc(sizeof(double) * 1000000);
double* b = (double*)malloc(sizeof(double) * 1000000);
double* c = (double*)malloc(sizeof(double) * 1000000);
for (int i = 0; i != 1000000; ++i)
{
double v;
scanf("%lf", &v);
a[i] = v;
b[i] = v;
c[i] = v;
}
gettimeofday(&t1, NULL);
mul_c(a, b);
gettimeofday(&t2, NULL);
time = 1000 * (t2.tv_sec - t1.tv_sec) + (t2.tv_usec - t1.tv_usec) / 1000;
printf("mul_c: %llu ms\n", time);
gettimeofday(&t1, NULL);
mul_asm(b, c);
gettimeofday(&t2, NULL);
time = 1000 * (t2.tv_sec - t1.tv_sec) + (t2.tv_usec - t1.tv_usec) / 1000;
printf("mul_asm: %llu ms\n\n", time);
for (int i = 0; i != 1000000; ++i)
{
printf("%lf\t\t\t%lf\n", a[i], b[i]);
}
return 0;
}
Я также пытался использовать все xmm регистрирует (0-7) и удаляет зависимости инструкций, чтобы улучшить параллельные вычисления:
void mul_asm(double* a, double* b)
{
asm volatile
(
".intel_syntax noprefix \n\t"
"xor rax, rax \n\t"
"0: \n\t"
"movupd xmm0, xmmword ptr [rdi+rax] \n\t"
"movupd xmm1, xmmword ptr [rdi+rax+16] \n\t"
"movupd xmm2, xmmword ptr [rdi+rax+32] \n\t"
"movupd xmm3, xmmword ptr [rdi+rax+48] \n\t"
"movupd xmm4, xmmword ptr [rdi+rax+64] \n\t"
"movupd xmm5, xmmword ptr [rdi+rax+80] \n\t"
"movupd xmm6, xmmword ptr [rdi+rax+96] \n\t"
"movupd xmm7, xmmword ptr [rdi+rax+112] \n\t"
"mulpd xmm0, xmmword ptr [rsi+rax] \n\t"
"mulpd xmm1, xmmword ptr [rsi+rax+16] \n\t"
"mulpd xmm2, xmmword ptr [rsi+rax+32] \n\t"
"mulpd xmm3, xmmword ptr [rsi+rax+48] \n\t"
"mulpd xmm4, xmmword ptr [rsi+rax+64] \n\t"
"mulpd xmm5, xmmword ptr [rsi+rax+80] \n\t"
"mulpd xmm6, xmmword ptr [rsi+rax+96] \n\t"
"mulpd xmm7, xmmword ptr [rsi+rax+112] \n\t"
"movupd xmmword ptr [rdi+rax], xmm0 \n\t"
"movupd xmmword ptr [rdi+rax+16], xmm1 \n\t"
"movupd xmmword ptr [rdi+rax+32], xmm2 \n\t"
"movupd xmmword ptr [rdi+rax+48], xmm3 \n\t"
"movupd xmmword ptr [rdi+rax+64], xmm4 \n\t"
"movupd xmmword ptr [rdi+rax+80], xmm5 \n\t"
"movupd xmmword ptr [rdi+rax+96], xmm6 \n\t"
"movupd xmmword ptr [rdi+rax+112], xmm7 \n\t"
"add rax, 128 \n\t"
"cmp rax, 8000000 \n\t"
"jne 0b \n\t"
".att_syntax noprefix \n\t"
:
: "D" (a), "S" (b)
: "memory", "cc"
);
}
Но он по-прежнему работает за 1 мс, с той же скоростью, что и обычная реализация C/C++.
ОБНОВЛЕНИЕ
Как следует из ответов / комментариев, я реализовал другой способ измерения времени выполнения:
#include <stdio.h>
#include <stdlib.h>
void mul_c(double* a, double* b)
{
for (int i = 0; i != 1000000; ++i)
{
a[i] = a[i] * b[i];
}
}
void mul_asm(double* a, double* b)
{
asm volatile
(
".intel_syntax noprefix \n\t"
"xor rax, rax \n\t"
"0: \n\t"
"movupd xmm0, xmmword ptr [rdi+rax] \n\t"
"mulpd xmm0, xmmword ptr [rsi+rax] \n\t"
"movupd xmmword ptr [rdi+rax], xmm0 \n\t"
"add rax, 16 \n\t"
"cmp rax, 8000000 \n\t"
"jne 0b \n\t"
".att_syntax noprefix \n\t"
:
: "D" (a), "S" (b)
: "memory", "cc"
);
}
void mul_asm2(double* a, double* b)
{
asm volatile
(
".intel_syntax noprefix \n\t"
"xor rax, rax \n\t"
"0: \n\t"
"movupd xmm0, xmmword ptr [rdi+rax] \n\t"
"movupd xmm1, xmmword ptr [rdi+rax+16] \n\t"
"movupd xmm2, xmmword ptr [rdi+rax+32] \n\t"
"movupd xmm3, xmmword ptr [rdi+rax+48] \n\t"
"movupd xmm4, xmmword ptr [rdi+rax+64] \n\t"
"movupd xmm5, xmmword ptr [rdi+rax+80] \n\t"
"movupd xmm6, xmmword ptr [rdi+rax+96] \n\t"
"movupd xmm7, xmmword ptr [rdi+rax+112] \n\t"
"mulpd xmm0, xmmword ptr [rsi+rax] \n\t"
"mulpd xmm1, xmmword ptr [rsi+rax+16] \n\t"
"mulpd xmm2, xmmword ptr [rsi+rax+32] \n\t"
"mulpd xmm3, xmmword ptr [rsi+rax+48] \n\t"
"mulpd xmm4, xmmword ptr [rsi+rax+64] \n\t"
"mulpd xmm5, xmmword ptr [rsi+rax+80] \n\t"
"mulpd xmm6, xmmword ptr [rsi+rax+96] \n\t"
"mulpd xmm7, xmmword ptr [rsi+rax+112] \n\t"
"movupd xmmword ptr [rdi+rax], xmm0 \n\t"
"movupd xmmword ptr [rdi+rax+16], xmm1 \n\t"
"movupd xmmword ptr [rdi+rax+32], xmm2 \n\t"
"movupd xmmword ptr [rdi+rax+48], xmm3 \n\t"
"movupd xmmword ptr [rdi+rax+64], xmm4 \n\t"
"movupd xmmword ptr [rdi+rax+80], xmm5 \n\t"
"movupd xmmword ptr [rdi+rax+96], xmm6 \n\t"
"movupd xmmword ptr [rdi+rax+112], xmm7 \n\t"
"add rax, 128 \n\t"
"cmp rax, 8000000 \n\t"
"jne 0b \n\t"
".att_syntax noprefix \n\t"
:
: "D" (a), "S" (b)
: "memory", "cc"
);
}
unsigned long timestamp()
{
unsigned long a;
asm volatile
(
".intel_syntax noprefix \n\t"
"xor rax, rax \n\t"
"xor rdx, rdx \n\t"
"RDTSCP \n\t"
"shl rdx, 32 \n\t"
"or rax, rdx \n\t"
".att_syntax noprefix \n\t"
: "=a" (a)
:
: "memory", "cc"
);
return a;
}
int main()
{
unsigned long t1;
unsigned long t2;
double* a;
double* b;
a = (double*)malloc(sizeof(double) * 1000000);
b = (double*)malloc(sizeof(double) * 1000000);
for (int i = 0; i != 1000000; ++i)
{
double v;
scanf("%lf", &v);
a[i] = v;
b[i] = v;
}
t1 = timestamp();
mul_c(a, b);
//mul_asm(a, b);
//mul_asm2(a, b);
t2 = timestamp();
printf("mul_c: %lu cycles\n\n", t2 - t1);
for (int i = 0; i != 1000000; ++i)
{
printf("%lf\t\t\t%lf\n", a[i], b[i]);
}
return 0;
}
Когда я запускаю программу с этим измерением, я получаю такой результат:
mul_c: ~2163971628 cycles
mul_asm: ~2532045184 cycles
mul_asm2: ~5230488 cycles <-- what???
Здесь стоит обратить внимание на две вещи, во-первых, количество циклов варьируется много, и я предполагаю, что это из-за операционной системы, позволяющей запускать другие процессы между ними. Есть ли способ предотвратить это или только считать циклы во время выполнения моей программы? Также, mul_asm2 производит идентичный выход по сравнению с двумя другими, но это намного быстрее, как?
Я попробовал программу Z boson в моей системе вместе со своими 2 реализациями и получил следующий результат:
> g++ -O2 -fopenmp main.cpp
> ./a.out
mul time 1.33, 18.08 GB/s
mul_SSE time 1.13, 21.24 GB/s
mul_SSE_NT time 1.51, 15.88 GB/s
mul_SSE_OMP time 0.79, 30.28 GB/s
mul_SSE_v2 time 1.12, 21.49 GB/s
mul_v2 time 1.26, 18.99 GB/s
mul_asm time 1.12, 21.50 GB/s
mul_asm2 time 1.09, 22.08 GB/s
3 ответа
В функции синхронизации, которую я использовал для предыдущих тестов, была серьезная ошибка. Это сильно недооценило пропускную способность без векторизации, а также других измерений. Кроме того, была другая проблема, которая заключалась в переоценке пропускной способности из-за COW в массиве, который был прочитан, но не записан. Наконец, максимальная пропускная способность, которую я использовал, была неправильной. Я обновил свой ответ с исправлениями и оставил старый ответ в конце этого ответа.
Ваша операция ограничена пропускной способностью памяти. Это означает, что процессор тратит большую часть своего времени на медленное чтение и запись в память. Отличное объяснение этому можно найти здесь: почему векторизация цикла не имеет улучшения производительности.
Однако я должен немного не согласиться с одним утверждением в этом ответе.
Таким образом, независимо от того, как это оптимизировано (векторизовано, развернуто и т. Д.), Оно не станет намного быстрее.
Фактически, векторизация , развертывание и множественные потоки могут значительно увеличить пропускную способность даже в операциях с ограниченной пропускной способностью памяти. Причина в том, что трудно получить максимальную пропускную способность памяти. Хорошее объяснение этому можно найти здесь: /questions/43644403/izmerenie-propusknoj-sposobnosti-pamyati-po-tochechnomu-proizvedeniyu-dvuh-massivov/43644415#43644415.
Остальная часть моего ответа покажет, как векторизация и несколько потоков могут приблизиться к максимальной пропускной способности памяти.
Моя тестовая система: Ubuntu 16.10, Skylake (i7-6700HQ@2.60GHz), 32 ГБ ОЗУ, двухканальный DDR4@2400 ГГц. Максимальная пропускная способность моей системы составляет 38,4 ГБ / с.
Из приведенного ниже кода я создаю следующие таблицы. Я установил номер потока, используя OMP_NUM_THREADS, например export OMP_NUM_THREADS=4, Эффективность bandwidth/max_bandwidth,
-O2 -march=native -fopenmp
Threads Efficiency
1 59.2%
2 76.6%
4 74.3%
8 70.7%
-O2 -march=native -fopenmp -funroll-loops
1 55.8%
2 76.5%
4 72.1%
8 72.2%
-O3 -march=native -fopenmp
1 63.9%
2 74.6%
4 63.9%
8 63.2%
-O3 -march=native -fopenmp -mprefer-avx128
1 67.8%
2 76.0%
4 63.9%
8 63.2%
-O3 -march=native -fopenmp -mprefer-avx128 -funroll-loops
1 68.8%
2 73.9%
4 69.0%
8 66.8%
После нескольких итераций бега из-за неопределенности измерений я сделал следующие выводы:
- Однопоточные скалярные операции получают более 50% пропускной способности.
- две многопоточные скалярные операции получают наибольшую пропускную способность.
- однопоточные векторные операции выполняются быстрее, чем однопоточные скалярные операции.
- однопоточные операции SSE быстрее, чем однопоточные операции AVX.
- Развертывание не помогает.
- Развертывание однопоточных операций выполняется медленнее, чем без развертывания.
- больше потоков, чем ядер (Hyper-Threading) обеспечивает меньшую пропускную способность.
Решение, которое обеспечивает лучшую пропускную способность, - это скалярные операции с двумя потоками.
Код, который я использовал для тестирования:
#include <stdlib.h>
#include <string.h>
#include <stdio.h>
#include <omp.h>
#define N 10000000
#define R 100
void mul(double *a, double *b) {
#pragma omp parallel for
for (int i = 0; i<N; i++) a[i] *= b[i];
}
int main() {
double maxbw = 2.4*2*8; // 2.4GHz * 2-channels * 64-bits * 1-byte/8-bits
double mem = 3*sizeof(double)*N*R*1E-9; // GB
double *a = (double*)malloc(sizeof *a * N);
double *b = (double*)malloc(sizeof *b * N);
//due to copy-on-write b must be initialized to get the correct bandwidth
//also, GCC will convert malloc + memset(0) to calloc so use memset(1)
memset(b, 1, sizeof *b * N);
double dtime = -omp_get_wtime();
for(int i=0; i<R; i++) mul(a,b);
dtime += omp_get_wtime();
printf("%.2f s, %.1f GB/s, %.1f%%\n", dtime, mem/dtime, 100*mem/dtime/maxbw);
free(a), free(b);
}
Старое решение с ошибкой синхронизации
Современное решение для встроенной сборки состоит в использовании встроенных функций. Есть еще случаи, когда требуется встроенная сборка, но это не один из них.
Одним из встроенных решений для вашего подхода к сборке в линейном режиме является просто:
void mul_SSE(double* a, double* b) {
for (int i = 0; i<N/2; i++)
_mm_store_pd(&a[2*i], _mm_mul_pd(_mm_load_pd(&a[2*i]),_mm_load_pd(&b[2*i])));
}
Позвольте мне определить тестовый код
#include <x86intrin.h>
#include <string.h>
#include <stdio.h>
#include <x86intrin.h>
#include <omp.h>
#define N 1000000
#define R 1000
typedef __attribute__(( aligned(32))) double aligned_double;
void (*fp)(aligned_double *a, aligned_double *b);
void mul(aligned_double* __restrict a, aligned_double* __restrict b) {
for (int i = 0; i<N; i++) a[i] *= b[i];
}
void mul_SSE(double* a, double* b) {
for (int i = 0; i<N/2; i++) _mm_store_pd(&a[2*i], _mm_mul_pd(_mm_load_pd(&a[2*i]),_mm_load_pd(&b[2*i])));
}
void mul_SSE_NT(double* a, double* b) {
for (int i = 0; i<N/2; i++) _mm_stream_pd(&a[2*i], _mm_mul_pd(_mm_load_pd(&a[2*i]),_mm_load_pd(&b[2*i])));
}
void mul_SSE_OMP(double* a, double* b) {
#pragma omp parallel for
for (int i = 0; i<N; i++) a[i] *= b[i];
}
void test(aligned_double *a, aligned_double *b, const char *name) {
double dtime;
const double mem = 3*sizeof(double)*N*R/1024/1024/1024;
const double maxbw = 34.1;
dtime = -omp_get_wtime();
for(int i=0; i<R; i++) fp(a,b);
dtime += omp_get_wtime();
printf("%s \t time %.2f s, %.1f GB/s, efficency %.1f%%\n", name, dtime, mem/dtime, 100*mem/dtime/maxbw);
}
int main() {
double *a = (double*)_mm_malloc(sizeof *a * N, 32);
double *b = (double*)_mm_malloc(sizeof *b * N, 32);
//b must be initialized to get the correct bandwidth!!!
memset(a, 1, sizeof *a * N);
memset(b, 1, sizeof *a * N);
fp = mul, test(a,b, "mul ");
fp = mul_SSE, test(a,b, "mul_SSE ");
fp = mul_SSE_NT, test(a,b, "mul_SSE_NT ");
fp = mul_SSE_OMP, test(a,b, "mul_SSE_OMP");
_mm_free(a), _mm_free(b);
}
Теперь первый тест
g++ -O2 -fopenmp test.cpp
./a.out
mul time 1.67 s, 13.1 GB/s, efficiency 38.5%
mul_SSE time 1.00 s, 21.9 GB/s, efficiency 64.3%
mul_SSE_NT time 1.05 s, 20.9 GB/s, efficiency 61.4%
mul_SSE_OMP time 0.74 s, 29.7 GB/s, efficiency 87.0%
Так с -O2 который не векторизует циклы, мы видим, что внутренняя версия SSE намного быстрее, чем простое решение C mul, efficiency = bandwith_measured/max_bandwidth где для моей системы максимальное значение составляет 34,1 ГБ / с.
Второй тест
g++ -O3 -fopenmp test.cpp
./a.out
mul time 1.05 s, 20.9 GB/s, efficiency 61.2%
mul_SSE time 0.99 s, 22.3 GB/s, efficiency 65.3%
mul_SSE_NT time 1.01 s, 21.7 GB/s, efficiency 63.7%
mul_SSE_OMP time 0.68 s, 32.5 GB/s, efficiency 95.2%
С -O3 Векторизация цикла, и встроенная функция не дает никаких преимуществ.
Третий тест
g++ -O3 -fopenmp -funroll-loops test.cpp
./a.out
mul time 0.85 s, 25.9 GB/s, efficency 76.1%
mul_SSE time 0.84 s, 26.2 GB/s, efficency 76.7%
mul_SSE_NT time 1.06 s, 20.8 GB/s, efficency 61.0%
mul_SSE_OMP time 0.76 s, 29.0 GB/s, efficency 85.0%
С -funroll-loops GCC развертывает циклы восемь раз, и мы видим значительное улучшение, за исключением решения для временного хранилища и не реальное преимущество для решения OpenMP.
Перед развертыванием петли сборка для mul Wiht -O3 является
xor eax, eax
.L2:
movupd xmm0, XMMWORD PTR [rsi+rax]
mulpd xmm0, XMMWORD PTR [rdi+rax]
movaps XMMWORD PTR [rdi+rax], xmm0
add rax, 16
cmp rax, 8000000
jne .L2
rep ret
С -O3 -funroll-loops сборка для mul является:
xor eax, eax
.L2:
movupd xmm0, XMMWORD PTR [rsi+rax]
movupd xmm1, XMMWORD PTR [rsi+16+rax]
mulpd xmm0, XMMWORD PTR [rdi+rax]
movupd xmm2, XMMWORD PTR [rsi+32+rax]
mulpd xmm1, XMMWORD PTR [rdi+16+rax]
movupd xmm3, XMMWORD PTR [rsi+48+rax]
mulpd xmm2, XMMWORD PTR [rdi+32+rax]
movupd xmm4, XMMWORD PTR [rsi+64+rax]
mulpd xmm3, XMMWORD PTR [rdi+48+rax]
movupd xmm5, XMMWORD PTR [rsi+80+rax]
mulpd xmm4, XMMWORD PTR [rdi+64+rax]
movupd xmm6, XMMWORD PTR [rsi+96+rax]
mulpd xmm5, XMMWORD PTR [rdi+80+rax]
movupd xmm7, XMMWORD PTR [rsi+112+rax]
mulpd xmm6, XMMWORD PTR [rdi+96+rax]
movaps XMMWORD PTR [rdi+rax], xmm0
mulpd xmm7, XMMWORD PTR [rdi+112+rax]
movaps XMMWORD PTR [rdi+16+rax], xmm1
movaps XMMWORD PTR [rdi+32+rax], xmm2
movaps XMMWORD PTR [rdi+48+rax], xmm3
movaps XMMWORD PTR [rdi+64+rax], xmm4
movaps XMMWORD PTR [rdi+80+rax], xmm5
movaps XMMWORD PTR [rdi+96+rax], xmm6
movaps XMMWORD PTR [rdi+112+rax], xmm7
sub rax, -128
cmp rax, 8000000
jne .L2
rep ret
Четвертый тест
g++ -O3 -fopenmp -mavx test.cpp
./a.out
mul time 0.87 s, 25.3 GB/s, efficiency 74.3%
mul_SSE time 0.88 s, 24.9 GB/s, efficiency 73.0%
mul_SSE_NT time 1.07 s, 20.6 GB/s, efficiency 60.5%
mul_SSE_OMP time 0.76 s, 29.0 GB/s, efficiency 85.2%
Теперь не встроенная функция является самой быстрой (исключая версию OpenMP).
Таким образом, в этом случае нет причин использовать встроенную или встроенную сборку, потому что мы можем добиться максимальной производительности с помощью соответствующих параметров компилятора (например, -O3, -funroll-loops, -mavx).
Тестовая система: Ubuntu 16.10, Skylake (i7-6700HQ@2.60GHz), 32 ГБ ОЗУ. Максимальная пропускная способность памяти (34,1 ГБ / с) https://ark.intel.com/products/88967/Intel-Core-i7-6700HQ-Processor-6M-Cache-up-to-3_50-GHz
Вот еще одно решение, которое стоит рассмотреть. cmp инструкция не нужна, если мы считаем от -N до нуля и обращаемся к массивам как N+i, GCC должен был это исправить давно. Это исключает одну инструкцию (хотя из-за макрооперации cmp и jmp часто считаются одной микрооперацией).
void mul_SSE_v2(double* a, double* b) {
for (ptrdiff_t i = -N; i<0; i+=2)
_mm_store_pd(&a[N + i], _mm_mul_pd(_mm_load_pd(&a[N + i]),_mm_load_pd(&b[N + i])));
Сборка с -O3
mul_SSE_v2(double*, double*):
mov rax, -1000000
.L9:
movapd xmm0, XMMWORD PTR [rdi+8000000+rax*8]
mulpd xmm0, XMMWORD PTR [rsi+8000000+rax*8]
movaps XMMWORD PTR [rdi+8000000+rax*8], xmm0
add rax, 2
jne .L9
rep ret
}
Эта оптимизация, возможно, будет полезна только для массивов, например, кеш L1, т.е. не для чтения из основной памяти.
Я наконец нашел способ получить простое C-решение, чтобы не генерировать cmp инструкция.
void mul_v2(aligned_double* __restrict a, aligned_double* __restrict b) {
for (int i = -N; i<0; i++) a[i] *= b[i];
}
А затем вызвать функцию из отдельного объектного файла, как это mul_v2(&a[N],&b[N]) так что это, пожалуй, лучшее решение. Однако если вы вызываете функцию из того же объектного файла (единицы перевода), который определен в GCC, генерируется cmp Инструкция снова.
Также,
void mul_v3(aligned_double* __restrict a, aligned_double* __restrict b) {
for (int i = -N; i<0; i++) a[N+i] *= b[N+i];
}
до сих пор генерирует cmp инструкция и генерирует ту же сборку, что и mul функция.
Функция mul_SSE_NT глупо Он использует невременные хранилища, которые полезны только при записи в память, но поскольку функция считывания и записи по одному и тому же адресу не временные хранилища, они не только бесполезны, но и дают худшие результаты.
Предыдущие версии этого ответа получали неверную пропускную способность. Причина была в том, что массивы не были инициализированы.
Ваш код asm действительно в порядке. Что не так, как вы это измеряете. Как я указал в комментариях, вы должны:
а) использовать больше итераций - 1 миллион - ничто для современного процессора
б) использовать HPT для измерения
c) использовать RDTSC или RDTSCP для подсчета реальных тактовых частот процессора
Кроме того, почему вы боитесь опций -O3? Не забудьте создать код для вашей платформы, поэтому используйте -march=native. Если ваш процессор поддерживает AVX или AVX2, компилятор получит возможность создавать еще лучший код.
Следующее - дайте компилятору несколько советов по псевдонимам и выравниванию, если вы знаете свой код.
Вот моя версия вашего mul_c - да, это специфично для GCC, но вы показали, что использовали GCC
void mul_c(double* restrict a, double* restrict b)
{
a = __builtin_assume_aligned (a, 16);
b = __builtin_assume_aligned (b, 16);
for (int i = 0; i != 1000000; ++i)
{
a[i] = a[i] * b[i];
}
}
Это произведет:
mul_c(double*, double*):
xor eax, eax
.L2:
movapd xmm0, XMMWORD PTR [rdi+rax]
mulpd xmm0, XMMWORD PTR [rsi+rax]
movaps XMMWORD PTR [rdi+rax], xmm0
add rax, 16
cmp rax, 8000000
jne .L2
rep ret
Если у вас AVX2 и выровнены данные на 32 байта, он станет
mul_c(double*, double*):
xor eax, eax
.L2:
vmovapd ymm0, YMMWORD PTR [rdi+rax]
vmulpd ymm0, ymm0, YMMWORD PTR [rsi+rax]
vmovapd YMMWORD PTR [rdi+rax], ymm0
add rax, 32
cmp rax, 8000000
jne .L2
vzeroupper
ret
Так что нет необходимости в asm, созданном вручную, если компилятор может сделать это для вас;)
Я хочу добавить еще одну точку зрения на проблему. Инструкции SIMD дают большой прирост производительности, если нет ограничений по объему памяти. Но в текущем примере слишком много операций загрузки и хранения памяти и слишком мало вычислений ЦП. Таким образом, процессор вовремя обрабатывает входящие данные без использования SIMD. Если вы используете данные другого типа (например, 32-разрядный с плавающей запятой) или более сложный алгоритм, пропускная способность памяти не будет ограничивать производительность процессора, и использование SIMD даст больше преимуществ.