Инициализация DSP-модуля

После того, как мы смогли корректно идентифицировать необходимые нам расширения процессора, возникает необходимость в инициализации всего DSP-модуля. Частично мы уже выполнили большую часть работы, написав код идентификации и поиска поддерживаемых расширений процессора. Теперь осталось открыть завесу над тайной, зачем всё это нам нужно. Давайте представим себе ситуацию, когда нам необходимо выполнить какую-то функцию, осуществляющую достаточно нетривиальные и интенсивные математические вычисления. Допустим, она будет иметь следующий прототип:

Код (C):

1.  
2. void example_function();
3.  

Очевидно, что для разных архитектур такая функция может иметь различные варианты реализации. Однако неизменным окажется то, что в любом случае потребуется написать нативную реализацию этой функции в случае, если возникнет необходимость скомпилировать код под «неоптимизированную» архитектуру. Поэтому первым делом пишем нативную реализацию этой функции в файле include/dsp/native/example.h:

Код (C):

1.  
2. #ifndef DSP_NATIVE_IMPL
3.     #error "This header should not be included directly"
4. #endif /* DSP_NATIVE_IMPL */
5.  
6. void example_function()
7. {
8.     printf("Native implementation of example function called\n");
9. }
10.  

При этом, реализацию этой функции организуем в отдельном пространстве имён dsp::native, а сам заголовочный файл включим внутри src/dsp/native.cpp:

Код (C):

1.  
2. #include <dsp.h>
3. #include <stdio.h>
4.  
5. #define DSP_NATIVE_IMPL
6.  
7. namespace dsp
8. {
9.     namespace native
10.     {
11.         #include <dsp/native/example.h>
12.     }
13. }
14.  
15. #undef DSP_NATIVE_IMPL
16.  

Разумеется, мы предполагаем, что у нас будет аналогичная оптимизированная функция для SSE и AVX, поэтому создаём файл include/dsp/x86/sse/example.h, где помещаем свою реализацию:

Код (C):

1.  
2. #ifndef X86_SSE_IMPL
3.     #error "This heades should not be included directly"
4. #endif /* X86_SSE_IMPL */
5.  
6. void example_function()
7. {
8.     printf("SSE-optimized implementation of example function called\n");
9. }
10.  

И аналогично в файле include/dsp/x86/avx/example.h:

Код (C):

1.  
2. #ifndef DSP_X86_AVX_IMPL
3.     #error "This header should not be included directly"
4. #endif /* DSP_X86_AVX_IMPL */
5.  
6. void example_function()
7. {
8.     printf("AVX-optimized implementation of example function called\n");
9. }
10.  

Аналогичным образом включаем эти файлы в src/dsp/sse.cpp:

Код (C):

1.  
2. #define X86_SSE_IMPL
3. namespace dsp
4. {
5.     namespace sse
6.     {
7.         #include <dsp/x86/sse/example.h>
8.     }
9. } /* namespace dsp::sse */
10.  
11. #undef X86_SSE_IMPL
12.  

И в src/dsp/avx.cpp:

Код (C):

1.  
2. #define DSP_X86_AVX_IMPL
3.  
4. namespace dsp
5. {
6.     namespace x86
7.     {
8.         #include <dsp/x86/avx/xcr.h>
9.     }
10. }
11.  
12. namespace dsp
13. {
14.     namespace avx
15.     {
16.         #include <dsp/x86/avx/example.h>
17.     }
18. }
19.  
20. #undef DSP_X86_AVX_IMPL
21.  

Однако теперь возникает закономерный вопрос – как приложение, использующее нашу библиотеку DSP, поймёт, какую функцию надо вызывать? Ведь в зависимости от архитектуры и конкретной модели процессора, одна реализация будет вести себя оптимально, другая – нет, а третья в принципе не может быть использована. Решить эту задачу можно, динамически формируя адрес оптимальной реализации функции. То есть, в области глобальных переменных мы объявим переменную-указатель на функцию, а вызывающий код будет читать этот адрес и передавать управление соответствующей реализации. Для этого в заголовочном файле include/dsp/dsp/example.h объявим соответствующий прототип переменной:

Код (C):

1.  
2. #ifndef DSP_H_IMPL
3.     #error "This header should not be included directly"
4. #endif /* DSP_H_IMPL */
5.  
6. /** This is an example function
7.  *
8.  */
9. extern void (* example_function)();
10.  

А саму переменную расположим, например, в src/dsp.cpp:

Код (C):

1.  
2. namespace dsp
3. {
4.     // Function prototypes
5.     void (* example_function)()     = NULL;
6. }
7.  

Однако указатель на функцию у нас изначально записан как NULL, и прежде чем иметь возможность вызвать оптимальную реализацию, необходимо этот указатель проинициализировать. Для этого тут же мы реализуем функцию инициализации DSP-модуля:

Код (C):

1.  
2. namespace dsp
3. {
4.     // Initialization
5.     void init()
6.     {
7.         native::init();
8.         #ifdef ARCH_X86
9.             x86::init();
10.         #endif /* ARCH_X86 */
11.     }
12. }
13.  

Как видно, нативная архитектура инициализируется в любом случае, а архитектура x86 – только в том случае, если исходный код компилируется под неё. Разумеется, файл src/dsp.cpp ничего не знает как о первой функции, так и о второй. Поэтому мы там же предварительно объявляем прототипы этих функций:

Код (C):

1.  
2. namespace dsp
3. {
4.     namespace native
5.     {
6.         void init();
7.     }
8.  
9. #ifdef ARCH_X86
10.     namespace x86
11.     {
12.         void init();
13.     }
14. #endif /* ARCH_X86 */
15. }
16.  

Теперь объявленные прототипы обязывают нас написать реализацию native::init в src/dsp/native.cpp:

Код (C):

1.  
2. namespace dsp
3. {
4.     namespace native
5.     {
6.         void init()
7.         {
8.             printf("Initializing native DSP functions\n");
9.             dsp::example_function       = native::example_function;
10.         }
11.     }
12. }
13.  

Ну а x86::init у нас уже частично реализована в src/dsp/x86.cpp, осталось только в самом конце вызвать функции инициализации расширений:

Код (C):

1.  
2. namespace dsp
3. {
4.     namespace x86
5.     {
6.         void init()
7.         {
8.             // Detect CPU options
9.             cpu_features_t f;
10.             detect_options(&f);
11.  
12.             // Extension initialization
13.             sse::init(&options);
14.             avx::init(&options);
15.         }
16.     }
17. }
18.  

Аналогичным образом чуть ранее объявляем прототипы для функций инициализации SSE и AVX:

Код (C):

1.  
2. namespace dsp
3. {
4.     namespace sse
5.     {
6.         void init(x86::cpu_features_t *options);
7.     }
8.  
9.     namespace avx
10.     {
11.         void init(x86::cpu_features_t *options);
12.     }
13. }
14.  

После этого в файле src/dsp/sse.cpp пишем реализацию sse::init:

Код (C):

1.  
2. namespace dsp
3. {
4.     namespace sse
5.     {
6.         using namespace x86;
7.  
8.         void init(cpu_features_t *options)
9.         {
10.             if (!(options->features & CPU_OPTION_SSE))
11.                 return;
12.  
13.             printf("Initializing SSE-optimized DSP functions\n");
14.  
15.             dsp::example_function       = sse::example_function;
16.         }
17.     }
18. }
19.  

И в файле src/dsp/avx.cpp – реализацию avx::init:

Код (C):

1.  
2. namespace dsp
3. {
4.     namespace avx
5.     {
6.         using namespace x86;
7.  
8.         void init(cpu_features_t *options)
9.         {
10.             if (!(options->features & CPU_OPTION_AVX))
11.                 return;
12.  
13.             printf("Initializing AVX-optimized DSP functions\n");
14.  
15.             dsp::example_function       = avx::example_function;
16.         }
17.     }
18. }
19.  

После всей проделанной работы нам необходимо где-то вызвать основную функцию инициализации – dsp::init. Открываем main.cpp и пишем в main:

Код (C):

1.  
2. #include <dsp.h>
3.  
4. int main()
5. {
6.     dsp::init();
7.  
8.     return 0;
9. }
10.  

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

Машинный контекст

Процессоры Intel в расширениях SSE и AVX поддерживают работу с числами с плавающей точкой стандарта IEEE 754. При этом, сами числа с плавающей точкой можно классифицировать на следующие виды:

• нули со знаком (Signed Zeros) – +0 и -0;
• обычные нормализованные числа;
• денормализованные числа;
• знаковые бесконечности – +Inf и -Inf;
• не числа (NaN, Not a Number) двух видов – SNaN и QNaN.

Особенность реализации вычислений над числами с плавающей точкой в архитектуре x86 заключается в том, что вычисления с участием денормализованных чисел приводит к большим скоростным регрессиям. Приведу пример из личного опыта: когда я тестировал только что разработанный мною эквалайзер, он замечательно и практически не загружая CPU справлялся со своей задачей при большом количестве одновременно включённых фильтров. Но стоило только мне убрать аудиосигнал, поступающий на вход, как эквалайзер начинал поедать ресурс CPU на все 100%. «Что же могло произойти?» – задумался я. На деле оказалось всё достаточно банально: отсутствие сигнала на входе эквалайзера равносильно подаче последовательности отсчётов с нулевой амплитудой, поэтому накапливаемые вещественные числа, хранящиеся в ячейках памяти цифровых фильтров, потихоньку начинают уменьшаться, и спустя весьма короткий промежуток времени нормализованные вещественные числа становятся денормализованными. Учитывая особенности реализация вещественной математики в процессорах семейства x86, имеем вполне закономерное последствие: производительность математических операций резко падает, а утилизация ресурсов центрального процессора подскакивает до 100%. Поэтому в нашем случае следует придерживаться одного жизненно необходимого правила: в критичных ко времени выполнения мультимедиа-приложениях следует избегать денормализованных чисел в расчётах.

В этом нам помогут специальные флаги FTZ и DAZ специального 32-разрядного регистра MXCSR, позволяющего контролировать и настраивать режим работы математики в SSE и AVX. Распишем значения его битов:

• бит 0 – IE (Invalid Operation Exception) – была выполнена неверная операция;
• бит 1 – DE (Denormal Exception) – была выполнена операция с денормализованным операндом;
• бит 2 – ZE (Divide-by-Zero Exception) – попытка деления на ноль;
• бит 3 – OE (Overflow Exception) – переполнение результата;
• бит 4 – UE (Underflow Exception) – антипереполнение результата;
• бит 5 – PE (Precision Exception) – потеря точности результата;
• бит 6 – DAZ (Denormals Are Zeros) – денормализованные операнды при расчёте приравниваются к нулям со знаком;
• бит 7 – IM (Invalid Operation Mask) – маскирование исключения IE;
• бит 8 – DM (Denormal Flag Mask) – маскирование исключения DE;
• бит 9 – ZM (Divide-by-Zero Mask) – маскирование исключения ZE;
• бит 10 – OM (Overflow Mask) – маскирование исключения OE;
• бит 11 – UM (Underflow Mask) – маскирование исключения UE;
• бит 12 – PM (Precision Mask) – маскирование исключения PE;
• биты 13 и 14 – RC (Rounding Control) – способ округления результата:
  • 00 – округление в сторону ближайшего, приоритет отдаётся чётному числу;
  • 01 – округление в сторону отрицательной бесконечности;
  • 10 – округление в сторону положительной бесконечности;
  • 11 – округление в сторону нуля;
• бит 15 – FTZ (Flush To Zero) – в случае, если установлен флаг OM и произошло антипереполнение результата, результат вычисления сбрасывается в ноль с сохранением знака, и устанавливаются флаги PE и UE.

Как видим, все биты регистра MXCSR управляют режимом работы вещественной арифметики, при этом использование битов FTZ и DAZ позволяет избегать денормализованных чисел за счёт частичной потери точности вычислений. Мы готовы этим пожертвовать, так как итоговая погрешность вычислений будет меньше пиковой амплитуды шума квантования 24-битного АЦП.

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

Согласно документации от Intel, по умолчанию для маски следует использовать значение 0xffbf, что, как раз, исключает оба бита FTZ и DAZ. Для того, чтобы получить актуальную маску для целевого процессора, необходимо воспользоваться специальной инструкцией сохранения контекста SSE, MMX и FPU – FXSAVE. Эта команда поддерживается только в том случае, если функция 1 инструкции CPUID возвращает в регистр ECX значение с установленным 26-м битом. Ранее, написав код детекции расширений процессора, мы учли это, и у нас есть отдельный бит индикации поддержки этой инструкции – CPU_OPTION_FXSAVE.

Теперь коротко о том, что делает FXSAVE: инструкция сохраняет регистры данных и управляющие регистры FPU, MMX и SSE в виде определённым образом структурированных данных в участке памяти размером 512 байт, при этом адрес участка памяти должен быть выровнен по границе 16 байт. Одним из элементов такой структуры, располагающимся по смещению 0x1c относительно начала и имеющим размер 4 байта, как раз, и является актуальная маска регистра MXCSR. Однако и тут кроются определённые подводные камни: дело в том, что далеко не все процессоры поддерживают запись маски MXCSR при выполнении FXSAVE, и значение поля следует дополнительно проверять. Если окажется, что оно равно нулю, то для MXCSR следует использовать маску по умолчанию.

Теперь мы знаем всё, что нам необходимо, и заводим дополнительный файл include/dsp/x86/sse/mxcsr.h, в котором реализуем функцию чтения маски MXCSR:

Код (C):

1.  
2. #ifndef CORE_X86_SSE_IMPL
3.     #error "This header should not be included directly"
4. #endif /* CORE_X86_SSE_IMPL */
5.  
6. #define MXCSR_IE                    (1 << 0)
7. #define MXCSR_DE                    (1 << 1)
8. #define MXCSR_ZE                    (1 << 2)
9. #define MXCSR_OE                    (1 << 3)
10. #define MXCSR_UE                    (1 << 4)
11. #define MXCSR_PE                    (1 << 5)
12. #define MXCSR_DAZ                   (1 << 6) /* Denormals are zeros flag */
13. #define MXCSR_IM                    (1 << 7)
14. #define MXCSR_DM                    (1 << 8)
15. #define MXCSR_ZM                    (1 << 9)
16. #define MXCSR_OM                    (1 << 10)
17. #define MXCSR_UM                    (1 << 11)
18. #define MXCSR_PM                    (1 << 12)
19. #define MXCSR_ALL_MASK              (MXCSR_IM | MXCSR_DM | MXCSR_ZM | MXCSR_OM | MXCSR_UM | MXCSR_PM)
20. #define MXCSR_RC_MASK               (3 << 13)
21. #define MXCSR_RC_NEAREST            (0 << 13)
22. #define MXCSR_RC_N_INF              (1 << 13)
23. #define MXCSR_RC_P_INF              (2 << 13)
24. #define MXCSR_RC_ZERO               (3 << 13)
25. #define MXCSR_FZ                    (1 << 15) /* Flush to zero flag */
26.  
27. #define MXCSR_DEFAULT               0xffbf
28.  
29. static uint32_t mxcsr_mask         = MXCSR_DEFAULT;
30.  
31. inline void init_mxcsr_mask()
32. {
33.     uint8_t fxsave[512] __lsp_aligned16;
34.     uint8_t *ptr        = fxsave;
35.  
36.     __asm__ __volatile__
37.     (
38.         // Clear FXSAVE structure
39.         __ASM_EMIT("xor     %%eax, %%eax")
40.         __ASM_EMIT("mov     $0x80, %%ecx")
41.         __ASM_EMIT("rep     stosl")
42.  
43.         // Issue fxsave
44.         __ASM_EMIT("fxsave  (%[fxsave])")
45.  
46.         // Get mask
47.         __ASM_EMIT("mov     0x1c(%[fxsave]), %%eax")
48.         __ASM_EMIT("test    %%eax, %%eax") // Old processors issue zero MXCSR mask
49.         __ASM_EMIT("jnz     1f")
50.         __ASM_EMIT("mov     %[mxcsr_dfl], %%eax")
51.         __ASM_EMIT("1:")
52.  
53.         // Store MXCSR mask
54.         __ASM_EMIT("mov     %%eax, %[mask]")
55.  
56.         : "+D" (ptr), [fxsave] "+S" (ptr), [mask] "+m" (mxcsr_mask)
57.         : [mxcsr_dfl] "i" (MXCSR_DEFAULT)
58.         : "cc", "memory",
59.           "%eax", "%ecx"
60.     );
61. }
62.  

Как видно из текста функции, в стеке выделяется структура fxsave, выровненная по границе 16 байт. После этого содержимое структуры заполняется нулями, после чего выполняется операция FXSAVE. Далее анализируется результат, сохранённый в структуре FXSAVE по смещению 0x1c, и если он равен нулю, то мы принимаем маску MXCSR равной MXCSR_DEFAULT (0xffbf). В противном случае мы сохраняем прочитанную из структуры FXSAVE маску в глобальной переменной mxcsr_mask.

Для того, чтобы иметь возможность устанавливать значение битов FTZ и DAZ в регистре MXCSR, нам следует написать функции чтения и записи регистра MXCSR:

Код (C):

1.  
2. inline uint32_t read_mxcsr()
3. {
4.     uint32_t result = 0;
5.  
6.     __asm__ __volatile__
7.     (
8.         __ASM_EMIT("stmxcsr %[result]")
9.         : [result] "+m" (result)
10.         :
11.         : "memory"
12.     );
13.  
14.     return result;
15. }
16.  
17. inline void write_mxcsr(uint32_t value)
18. {
19.     __asm__ __volatile__
20.     (
21.         // Clear FXSAVE structure
22.         __ASM_EMIT("and         %[mask], %[value]")
23.         __ASM_EMIT("ldmxcsr     %[value]")
24.         : [value] "+m" (value)
25.         : [mask] "r" (mxcsr_mask)
26.         : "cc", "memory"
27.     );
28. }
29.  

Отлично, теперь мы можем подключить функции управления регистром в файл src/dsp/sse.cpp:

Код (C):

1.  
2. #define X86_SSE_IMPL
3. namespace dsp
4. {
5.     namespace sse
6.     {
7.         #include <dsp/x86/sse/mxcsr.h>
8.     }
9. } /* namespace dsp::sse */
10.  
11. #undef X86_SSE_IMPL
12.  

Разумеется, хорошим тоном будет, если в случае необходимости выполнения математических расчётов мы будем по возможности устанавливать флаги DAZ и FTZ в регистре MXCSR, а после выполнения расчётов – восстанавливать исходное состояние регистра MXCSR. Для этого введём понятие машинного контекста в файле include/dsp.h:

Код (C):

1.  
2. namespace dsp
3. {
4.     #pragma pack(push, 1)
5.     typedef struct context_t
6.     {
7.         uint32_t        top;
8.         uint32_t        data[15];
9.     } context_t;
10.     #pragma pack(pop)
11. }
12.  

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

Код (C):

1.  
2. namespace dsp
3. {
4.     // Start and finish types
5.     typedef void (* start_t)(context_t *ctx);
6.     typedef void (* finish_t)(context_t *ctx);
7.  
8.     /** Start DSP processing, save machine context
9.      *
10.      * @param ctx structure to save context
11.      */
12.     extern void (* start)(context_t *ctx);
13.  
14.     /** Finish DSP processing, restore machine context
15.      *
16.      * @param ctx structure to restore context
17.      */
18.     extern void (* finish)(context_t *ctx);
19.  
20. }
21.  

Как видно, start и stop – указатели на функции, которые должны быть инициализированы в ходе выполнения функции dsp::init. Делаем реализацию этих функций для нативной архитектуры:

Код (C):

1.  
2. namespace dsp
3. {
4.     namespace native
5.     {
6.         void start(context_t *ctx)
7.         {
8.             ctx->top        = 0;
9.         }
10.  
11.         void finish(context_t *ctx)
12.         {
13.         }
14.  
15.         void init()
16.         {
17.             printf("Initializing native DSP functions\n");
18.             dsp::start                  = native::start;
19.             dsp::finish                 = native::finish;
20.             dsp::example_function       = native::example_function;
21.         }
22.     }
23. }
24.  

Как видно из реализации, функция start инициализирует вершину стека, а finish ничего не делает ввиду того, что восстанавливать, особо, ничего и не надо.

В случае использования SSE нам потребуется сохранить адреса этих функций и переназначить на функции, которые будут изменять значение регистра MXCSR. Сами же переназначенные функции должны будут вызвать сохранённые указатели в прологе и эпилоге соответственно:

Код (C):

1.  
2. namespace dsp
3. {
4.     namespace sse
5.     {
6.         using namespace x86;
7.  
8.         static dsp::start_t     dsp_start       = NULL;
9.         static dsp::finish_t    dsp_finish      = NULL;
10.  
11.         static void start(context_t *ctx)
12.         {
13.             dsp_start(ctx);
14.             uint32_t    mxcsr       = read_mxcsr();
15.             ctx->data[ctx->top++]   = mxcsr;
16.             write_mxcsr(mxcsr | MXCSR_FZ | MXCSR_DAZ);
17.         }
18.  
19.         static void finish(context_t *ctx)
20.         {
21.             write_mxcsr(ctx->data[--ctx->top]);
22.             dsp_finish(ctx);
23.         }
24.  
25.         void init(cpu_features_t *options)
26.         {
27.             if (!(options->features & CPU_OPTION_SSE))
28.                 return;
29.  
30.             printf("Initializing SSE-optimized DSP functions\n");
31.             dsp_start                   = dsp::start;
32.             dsp_finish                  = dsp::finish;
33.             dsp::start                  = sse::start;
34.             dsp::finish                 = sse::finish;
35.  
36.             dsp::example_function       = sse::example_function;
37.         }
38.     }
39. }
40.  

Теперь, для того, чтобы вызвать нашу example_function, мы должны следовать определённой конвенции: вызывать функцию dsp::start вначале и dsp::finish в конце:

Код (C):

1.  
2. void do_math()
3. {
4.     dsp::context_t ctx;
5.     dsp::start(&ctx);
6.  
7.     dsp::example_function();
8.  
9.     dsp::finish(&ctx);
10. }
11.  
12. int main()
13. {
14.     dsp::init();
15.  
16.     do_math();
17.  
18.     return 0;
19. }
20.  

Сборка проекта

Осталось собрать проект. В данном случае пишем обычный Makefile, который для сборки бинарников использует отдельный каталог. Отдельное внимание уделим мейкфайлу в каталоге src/dsp:

Код (Text):

1.  
2. FILE                    = $(@:$(OBJDIR)/%.o=%.cpp)
3.  
4. ALL_IMPL                = $(OBJDIR)/dsp.o
5. NATIVE_IMPL             = $(OBJDIR)/native.o
6. X86_IMPL                = $(OBJDIR)/x86.o
7. SSE_IMPL                = $(OBJDIR)/sse.o
8. AVX_IMPL                = $(OBJDIR)/avx.o
9.  
10. SSE_INSTR_FLAGS         = -mmmx -m3dnow -msse
11. AVX_INSTR_FLAGS         = -mavx
12.  
13. # Form set of object to compile and link
14. LINK_OBJECTS            = $(NATIVE_IMPL)
15. ifeq ($(CPU_ARCH), i586)
16. LINK_OBJECTS           += $(X86_IMPL) $(SSE_IMPL) $(AVX_IMPL)
17. endif
18. ifeq ($(CPU_ARCH), x86_64)
19. LINK_OBJECTS           += $(X86_IMPL) $(SSE_IMPL) $(AVX_IMPL)
20. endif
21.  
22. .PHONY: all
23.  
24. all: $(ALL_IMPL)
25.  
26. $(ALL_IMPL): $(LINK_OBJECTS)
27.     @echo "  $(LD) $(notdir $(ALL_IMPL))"
28.     @$(LD) -r $(LDFLAGS) -o $(ALL_IMPL) $(LINK_OBJECTS)
29.  
30. $(NATIVE_IMPL) $(X86_IMPL):
31.     @echo "  $(CC) $(FILE)"
32.     @$(CC) -c $(CPPFLAGS) $(CFLAGS) $(INCLUDE) $(FILE) -o $(@)
33.  
34. $(SSE_IMPL):
35.     @echo "  $(CC) $(FILE)"
36.     @$(CC) -c $(CPPFLAGS) $(CFLAGS) $(SSE_INSTR_FLAGS) $(INCLUDE) $(FILE) -o $(@)
37.  
38. $(AVX_IMPL):
39.     @echo "  $(CC) $(FILE)"
40.     @$(CC) -c $(CPPFLAGS) $(CFLAGS) $(AVX_INSTR_FLAGS) $(INCLUDE) $(FILE) -o $(@)
41.  

Как видно, отдельные файлы компилируются с отдельными ключами, включающими те или иные расширенные наборы инструкций. Если мы в GCC попытаемся скомпилировать SSE-код с ключом -mavx, то компилятор автоматически заменит все SSE-инструкции на их AVX-аналоги, что, в итоге, может привести к некорректно работающему коду ввиду того, что детали реализации отдельных инструкций в AVX отличаются от деталей реализации аналогичных им SSE-инструкций. Поэтому код, использующий SSE должен компилироваться с ключами для SSE, а код, использующий AVX – с ключами для AVX.

В итоге на выходе мы получаем каталог .build, в котором располагается целевой файл dsp-utils. Запустим его:

Код (Text):

1.  
2. .build/dsp-utils
3. Initializing native DSP functions
4. Initializing SSE-optimized DSP functions
5. Initializing AVX-optimized DSP functions
6. AVX-optimized implementation of example function called
7.  

Как видно, наш код инициализации обнаружил расширение AVX на целевой архитектуре и выбрал в качестве оптимального решения AVX-реализацию example_function.

Давайте посмотрим, как поведёт себя код на процессоре, не поддерживающем AVX:

Код (Text):

1.  
2. .build/dsp-utils
3. Initializing native DSP functions
4. Initializing SSE-optimized DSP functions
5. SSE-optimized implementation of example function called
6.  

Как видно, мы получили вполне себе ожидаемый результат.

К сожалению, компьютер, не поддерживающий SSE, сейчас найти трудно, так что будем полагать, что он корректно работает и в случае с отсутствием поддержки SSE.

Заключение

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

Я думал, как лучше всего организовать хранение примеров к данному циклу статей и пришёл к выводу, что лучше всего с этим справится система контроля версий. Поэтому я завёл отдельный Git-репозиторий для хранения исходных кодов:
https://github.com/sadko4u/simd-utils

При этом, чтобы не возникало путаницы, я решил тэговать определённый коммит к моменту выхода определённой статьи. Поэтому актуальный для данной статьи исходный код доступен по следующему адресу:
https://github.com/sadko4u/simd-utils/tree/SIMD-003

Если же вы вытащили master-ветку, то переключиться на соответствующий тег не составит большой проблемы:

Код (Text):

1.  
2. git checkout tags/SIMD-003
3.