Элементы оптимизации прикладных программ для Intel Xeon Phi. Intel C/C++ Compiler

Использование дополнительных параметров offload директивы применительно к предыдущему примеру приведено ниже:

__attribute__((target(mic))) void func(float* a, 
    float* b, int count, float c, float d)
{
    #pragma omp parallel for 
    for (int i = 0; i < count; ++i)
    {
        a[i] = b[i]*c + d;
    }
}

int main()
{
    const int count = 100;
    float a[count], b[count], c, d;
    …
    #pragma offload target(mic) in(b) out(a)
        func(a, b, count, c, d);
    …
}

Теперь выполняется обмен только теми данными, которые нужны. Длина массивов здесь не указывается, т.к. массивы статические. Переменные count, c и d будут скопированы на сопроцессор автоматически.

Работа с динамическими массивами чуть сложнее. В частности, необходимо явно указывать размер данных при копировании:

#define ALLOC alloc_if(1) free_if(0)
#define FREE alloc_if(0) free_if(1)
#define REUSE alloc_if(0) free_if(0)

void f()
{
    int *p = (int *)malloc(100*sizeof(int));
    // Memory is allocated for p, 
    // data is sent from CPU and retained
    #pragma offload target(mic:0) in(p[0:100] : ALLOC)
    { p[6] = 66; }
    …
    // Memory for p reused from previous offload 
    // and retained once again
    // Fresh data is sent into the memory
    #pragma offload target(mic:0) in(p[0:100] : REUSE)
    { p[6] = 66; }
    …
    // Memory for p reused from previous offload,
    // freed after this offload.
    // Final data is pulled from coprocessor to CPU
    #pragma offload target(mic:0) out(p[0:100] : FREE)
    { p[7] = 77; }
    …
}

В данном примере демонстрируется также возможность хранения данных на сопроцессоре между вызовами кода на нем.

Рассмотрим этот пример подробней.

Во-первых, мы динамически выделяем память на центральном процессоре. К этому моменту данные существуют только в оперативной памяти CPU.

Первое обращение к сопроцессору приводит к выделению памяти на нем и копированию туда данных из оперативной памяти. После чего на MIC изменяется элемент p[6] . После окончания этого шага имеем почти идентичные копии массива p[] на CPU и сопроцессоре. Отличается один элемент p[6] , который равен 66 на сопроцессоре. По окончании шага данные на сопроцессоре не удаляются.

При втором обращении к сопроцессору выделения/удаления памяти не происходит, однако происходит копирование данных из оперативной памяти в память сопроцессора, т.е. выполняется обновление массива p[].

На заключительном шаге данные с сопроцессора копируются в оперативную память, после чего память на MIC освобождается. Обратите внимание, что элемент p[7] на CPU будет равен 77 по завершении операции.

Рассмотрим еще один пример работы с динамической памятью. Но на этот раз будем использовать данные, размещенные только на сопроцессоре:

void f()
{
    int *p;
    … 
    // The nocopy clause ensures CPU values pointed to by p
    // are not transferred to coprocessor
    #pragma offload target(mic:0) nocopy(p)
    {
        // Allocate dynamic memory for p on coprocessor
        p = (int *)malloc(100);
        p[0] = 77;
        …
    }
    ..
    // The nocopy clause ensures p is not altered 
    // by the offload process
    #pragma offload target(mic:0) nocopy(p)
    {
        // Reuse dynamic memory pointed to by p
        … = p[0]; // Will be 77
    }

}

Обратите внимание, что при использовании модификатора nocopy() автоматического выделения, а значит и удаления памяти не происходит. Это позволяет нам выделять память только в рамках сопроцессора и использовать ее для хранения данных между вызовами кода для Intel MIC.

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

Приведенный ниже пример показывает, как организовать одновременную работу процессора и сопроцессора:

double __attribute__((target(mic)))
myworkload(double input)
{
    // do something useful here
    return result;
}

int main(void)
{
    //…. Initialize variables
    #pragma omp parallel sections
    {
        #pragma omp section
        {
            #pragma offload target(mic) 
                result1= myworkload(input1);
        }  
        #pragma omp section
        {
            result2= myworkload(input2);
        }
    }
}

Для создания вспомогательных CPU-потоков не обязательно использовать возможности технологии OpenMP. Подойдет любое доступное вам средство, начиная от потоков операционной системы и заканчивая специальными библиотеками (например, TBB).

Следующий пример демонстрирует применение асинхронного offload режима для организации схемы двойной буферизации данных. Метод двойной буферизации применяется для того, чтобы уменьшить время передачи данных за счет одновременного выполнения передач и полезных расчетов (рис 5.2).

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

В коде с использованием асинхронных возможностей offload режима это выглядит следующим образом:

int main(int argc, char* argv[])
{
    // Allocate  initialize in1, res1,
    // in2, res2 on host
    #pragma offload_transfer target(mic:0) in(cnt)\
        nocopy(in1, res1, in2, res2 : length(cnt) \
        alloc_if(1) free_if(0))

    do_async_in();

    // Free MIC memory
    #pragma offload_transfer target(mic:0) \
        nocopy(in1, res1, in2, res2 : length(cnt) \
        alloc_if(0) free_if(1))

    return 0;
}

void do_async_in() 
{
    float lsum;
    int i;
    lsum = 0.0f;
    
    #pragma offload_transfer target(mic:0) \
        in(in1 : length(cnt) \
        alloc_if(0) free_if(0)) signal(in1)
    for (i = 0; i < iter; i++) 
    {
        if (i % 2 == 0) 
        {
            #pragma offload_transfer target(mic:0) \
                if(i !=iter - 1) in(in2 : length(cnt) \
                alloc_if(0) free_if(0)) signal(in2)
            
            #pragma offload target(mic:0) nocopy(in1) \
                wait(in1) out(res1 : length(cnt) \
                alloc_if(0) free_if(0))
            {
                compute(in1, res1);
            }
            
            lsum = lsum + sum_array(res1);
        } 
        else 
        {
            #pragma offload_transfer target(mic:0) \
                if(i != iter - 1) in(in1 : length(cnt) \
                alloc_if(0) free_if(0))  signal(in1)

            #pragma offload target(mic:0) nocopy(in2) \
                wait(in2) out(res2 : length(cnt) \
                alloc_if(0) free_if(0))
            {
                compute(in2, res2);
            }
            
            lsum = lsum + sum_array(res2);
        }
    } // for
} // do_async_in()

В завершении раздела о явной схеме работы с памятью рассмотрим пример работы со сложными структурами данных и дадим метод их передачи на сопроцессор:

typedef struct {
    int m1;
    char *m2;
} nbwcs;

void sample11()
{
    nbwcs struct1;
    struct1.m1 = 10;
    struct1.m2 = malloc(11);
    int m1;
    char *m2;
    
    // Disassemble the struct for transfer to target
    m1 = struct1.m1;
    m2 = struct1.m2;
    
    #pragma offload target(mic) inout(m1) inout(m2 : length(11)
    {
        nbwcs struct2;
        // Reassemble the struct on the target
        struct2.m1 = m1;
        struct2.m2 = m2;
        ...
    }
    ...
}

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

Неявная схема работы с памятью в режиме offload

Вторым способом работы с сопроцессором является неявный режим, основанный на использовании ключевых слов – расширений языков C/C++ (не поддерживается языком Fortran). Основная идея неявной схемы состоит в использовании разделяемой между CPU и MIC памяти в рамках единого виртуального адресного пространства (рис 5.3).

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

Для выделения участка разделяемой памяти динамически необходимо воспользоваться функциями:

Удаление памяти выполняется соответственно с помощью функций:

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

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

Заметим, что неявная схема является частью расширения Intel Cilk Plus.

Ключевое слово _Cilk_shared используется для размещения переменных и функций в разделяемой памяти. Примеры его использования приведены в таблице 5.4.

 int _Cilk_shared f(int x)
{ return x + 1; } 

 _Cilk_shared int x = 0; 

 static _Cilk_shared int x; 

 class _Cilk_shared x {…}; 

 int _Cilk_shared *p; 

 int* _Cilk_shared p; 

 #pragma offload_attribute(push, target(mic))
…
#pragma offload_attribute(pop) 

Выполнение кода на сопроцессоре обеспечивается ключевым словом _Cilk_offload (таблица 5.5).

 x = _Cilk_offload func(y); 

 x = _Cilk_offload_to(card_num) func(y); 

 x = _Cilk_spawn _Cilk_offload func(y); 

 _Cilk_offload _Cilk_for(i=0; i<N; ++i)
{ a[i] = b[i] + c[i]; } 

Использование неявной схемы продемонстрируем на примере функции для вычисления числа Пи (рис 5.4).

В данном примере используется глобальная переменная pi, размещаемая в разделяемой памяти CPU и сопроцессора. Используется функция compute_pi(), объявленная как общая для CPU и MIC. При вызове указывается необходимость ее запуска на сопроцессоре, если это возможно.

Сравнение явной и неявной схемы работы с памятью

В таблице 6 приводятся основные моменты, на которые стоит обратить внимание при выборе той или иной схемы работы с памятью. Вкратце, неявная схема более удобна в использовании, а явная предоставляет больше возможностей для контроля и оптимизации.

Для обеих схем компилятор генерирует два типа бинарных кодов процессорную и сопроцессорную версию. Версия для CPU содержит все переменные и функции независимо от того, отмечены ли они offload директивами (ключевыми словами), или нет. Сопроцессорная версия содержит только функции и переменные, отмеченные как offload. Обе версии объединены в один исполняемый файл.

Векторизация

В этом разделе рассматриваются возможности векторизации приложений с использованием Intel C/C++ Compiler. Информация, представленная в данном разделе, актуальна для всех моделей программирования на сопроцессоре, а также применима при программировании на обычном CPU.

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

Рассмотрим простой пример:

float *restrict A, *B, *C;
for (int i = 0; i < n; ++i) 
{
    A[i] = B[i] + C[i];
}

При исполнении такого кода в скалярном виде процессор будет выполнять одно сложение за такт и потратит на этот участок n тактов. В то же время современный процессор с поддержкой SSE может за такт выполнить 4 сложения, с поддержкой AVX – 8, а сопроцессор Intel Xeon Phi – 16 (рис 5.5). А это означает, что если такой код будет скомпилирован с использованием векторных инструкций процессора, то он может быть выполнен в несколько раз быстрее.

Каким образом можно сделать свой код векторным, используя компилятор компании Intel? Существует несколько вариантов:

В данной лекции будут рассмотрены только две первые возможности - автоматическая векторизация и применение возможностей Intel Cilk Plus.

Автоматическая векторизация

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

Рассмотрим процесс векторизации более подробно. В качестве примера возьмем следующий код:

for(i=0;i<*p;i++)
{
    A[i] = B[i]*C[i];
    sum = sum + A[i];
}

Перед тем, как выполнить автоматическую векторизацию кода, компилятор пытается проверить выполнение следующих условий:

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

Например, для компилятора часто сложным является ответ на вопрос о том, что массив A[] не перекрываются с массивами B[] и C[]. Для того чтобы отразить это в синтаксисе языка, можно объявить указатель A с ключевым словом restrict:

float* restrict A;

Такое объявление говорит компилятору о том, что массив A[] не перекрывается с другими массивами.

Вопрос определения инвариантов цикла для компилятора тоже не тривиален. По умолчанию если компилятор не может принять решения о том, является ли та или иная переменная является инвариантом, он считает ее не инвариантом, тем самым обеспечивая корректность кода. Однако и векторизацию такого цикла компилятор выполнить не может. Для того чтобы сказать компилятору об отсутствии зависимостей в цикле, используется директива #pragma ivdep перед телом цикла:

#pragma ivdep
for(i=0; i < *p; i++)
{
    A[i] = B[i]*C[i];
    sum = sum + A[i];
}

Директива ivdep дает команду компилятору игнорировать недоказанные зависимости. Однако если компилятор нашел и доказал зависимость, то векторизации цикла даже с этой директивой не произойдет.

Отметим, что данная директива поддерживаются и языком Fortran, где имеет вид !dir$ ivdep.

Использование директивы SIMD

Векторизация с помощью директивы #pragma simd дополняет автоматическую векторизацию так же, как распараллеливание с помощью #pragma omp дополняет автораспараллеливание. По сути, директивы simd и omp являются прямыми аналогами, позволяя выполнить векторизацию и распараллеливание вручную. При этом корректность работы программы в обоих случая не гарантируется и должна обеспечиваться разработчиком.

Пример использования директивы simd приведен ниже [5.3]:

void add_floats(float *a, float *b, float *c, float *d, float *e, int n){
    int i;
    
    #pragma simd
    for (i=0; i<n; i++){
        a[i] = a[i] + b[i] + c[i] + d[i] + e[i];
    }
}

Как и #pragma omp, директива simd может содержать дополнительные параметры, посредством которых можно сообщить компилятору о том, как корректно и эффективно векторизовать данный участок кода. Полное описание директивы содержится в соответствующем разделе документации по компилятору [5.4].

Рассмотрим основные параметры simd директивы:

#pragma simd vectorlength(4)
for (i = 0; i < n; i++) {
    a[i] = a[i] + b[i] + c[i];
}

#pragma simd linear(k:j)
for (i = 0; i < n; i += step) {
    k += j;
    a[i] = a[i] + b[n - k + 1];
}

int x = 0;
#pragma simd reduction(+:x)
for (i = 0; i < n; ++i)
    x = x + A[i];

double pi(int count) 
{
    int i;
    double pi = 0.0;
    double t;
    #pragma simd private(t) reduction(+:pi)
    for (i=0; i<count; i++) {
        t = (double)((i+0.5)/count);
        pi += 4.0/(1.0+t*t);
    }
    pi /= count;
    return pi;
}

Отметим, что на настоящий момент в компиляторах Intel активно улучшается поддержка процесса векторизации с использованием директивы simd, повышается область ее применения и добавляются новые возможности.

Технология Array Notation в Intel Cilk Plus

Другой способ векторизации кода – использование расширения Intel Cilk Plus для языков программирования C/C++ в части векторизации кода. Для нашего примера лучше всего подойдет применение специальной технологии Array Notation, с использованием которой код можно переписать так (обратите внимание на отсутствие цикла):

A[:] = B[:] + C[:];

Использование специальной технологии Array Notation является мощным инструментом для векторизации и значительно более сложных участков кода. Остановимся на ней подробнее.

Во-первых, обсудим синтаксис этой нотации:

Многомерные массивы также поддерживаются.

Во-вторых, приведем список действий, которые можно применять к массивам в этой нотации:

d[:] = a[:] + (b[:]*c[:]);

b[:] = func(a[:]);

sum = __sec_reduce_add(a[:]);

if (mask[:])
{
    a[:] = b[:];
}

c[:] = a[b[:]]; //gather
a[b[:]] = c[:]; //scatter

b[:] = __sec_shift_right(a[:], shift_val, fill_val);
b[:] = __sec_shift_left(a[:], shift_val, fill_val);
b[:] = __sec_rotate_right(a[:], rotate_val);
b[:] = __sec_rotate_left(a[:], rotate_val);

И наконец, опишем принципы работы этого механизма:

a[0:5] = b[0:6]; // No. Size mismatch.
a[0:5][0:4] = b[0:5]; // No. Rank mismatch.
a[0:5] = b[0:5][0:5]; // No. No 2D->1D reduction.
a[0:4] = 5; // OK. 4 elements of A filled w/ 5.
a[0:4] = b[i]; // OK. Fill with scalar b[i].
a[10][0:4] = b[1:4]; // OK. Both are 1D sections.
b[i] = a[0:4]; // No. Use reduction intrinsic.

Приведем еще один пример использования специальной технологии Array Notation. Выполняется скалярное произведение векторов. Скалярный код выглядит так:

float dot_product(unsigned int size, float *A, float *B)
{
    int i;
    float dp=0.0f;
    for (i=0; i<size; i++) 
    {
        dp += A[i] * B[i];
    }
    
    return dp;
}

С использованием Intel Cilk Plus его можно написать так:

float dot_product(unsigned int size, float A[size], float B[size])
{
    return __sec_reduce_add(A[:] * B[:]);
}

В обеих версиях используются статические массивы размера size. При этом вторая версия кода будет работать с использованием векторных расширений процессора или сопроцессора.

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

#include <cilk\cilk.h>

void saxpy_vec(int m, float a, float x[m], float y[m]) 
{
    y[0:m] += a*x[0:m];
}

void main(void) 
{
    int n = 2048;
    const int m = 256;
    float* a = new float[n];
    float* b = new float[n];
    cilk_for (int i = 0; i < n; i += m) 
    {
        saxpy_vec(m, 2.0, a[i], b[i]);
    }
}

Элементарные функции в Intel Cilk Plus

Использование механизма элементарных функций в Intel Cilk Plus является еще одной возможностью векторизации вашего кода.

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

__declspec(vector) float foo(float *B, float *C, int i)
{
    return B[i] + C[i];
}
…
for(i=0; i<N; i++)
{
    A[i] = foo(B, C, i);
}
…

Основная идея здесь состоит в выделении векторизуемой операции в отдельную функцию и вызове этой функции в цикле.

Рассмотрим этот механизм подробнее.

Элементарная функция – это функция, выполняющая вычисления над скалярными элементами данных. Строка __declspec(vector) указывает на необходимость векторизации кода этой функции:

__declspec(vector) float foo(float a, float b,
    float c, float d) 
{
    return a * b + c * d;
}

Вызывать элементарную функцию можно следующим образом:

for (i=0; i<n; i++)
{
    A[i] = foo(B[i], C[i], D[i], E[i]);
}

A[:] = foo(B[:], C[:], D[:], E[:]);

__declspec(vector) float bar(float a, float b,
    float c, float d)
{
    return sinf(foo(a,b,c,d));
}

e = foo(a, b, c, d);

При этом во всех случаях кроме последнего будет получен векторизованный код.

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

Для обеспечения векторизации кода с использованием элементарных функций программист должен написать обычную C/C++ скалярную функцию, добавить к ней описание __declspec(vector) и вызвать ее в параллельном контексте (одним из трех способов, описанных выше).

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

Заметим, что для правильной работы механизма элементарных функций следует соблюдать следующие правила:

Подробнее об использовании конструкций Intel Cilk Plus можно узнать по ссылкам [5.14],[5.15].

Методы эффективной векторизации кода с использованием Intel C/C++ Compiler

Использование отчетов компилятора

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

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

Опция компилятора, отвечающая за отчеты по векторизации, имеет вид -vec-report[n] для Linux и /Qvec-report[n] для Windows версии компилятора.

Здесь n может принимать значения от 0 до 7. По умолчанию используется значение 0, что означает отсутствие каких-либо сообщений. До недавнего времени максимум информации давала опция с n равным 3. А именно, с ее помощью можно было понять, какие циклы успешно векторизовались, а какие нет и почему.

Недавнее появление режима 6 позволяет существенно расширить эти данные за счет более конкретной информации о причинах неудачи компилятора (режим 6 появился в Intel Composer 13) [5.7]. Теперь компилятор говорит не только о том, почему у него не получилось векторизовать цикл, но и о том, что нужно сделать пользователю, чтобы исправить ситуацию.

Например, при использовании опции –vec-report3 вы можете получить сообщение вида:

loop was not vectorized: unsupported data type

А использование –vec-report6 даст уже:

vectorization support: type TTT is not supported for operation OOO

Т.е. разработчик может попытаться изменить тип данных в указанной операции или саму операцию для автоматической векторизации.

Опция –vec-report7 появилась уже в Intel Composer версии 14. Она позволяет получить дополнительную информацию о векторизуемом коде, например, ожидаемое ускорение от векторизации, используемый здесь шаблон доступа к памяти и др.

Выравнивание данных

При работе с векторными регистрами следует обращать внимание на эффективность чтения/записи данных в эти регистры. В силу особенностей архитектуры процессора, операции чтения и записи наиболее эффективно работают с данными, которые являются выровненными по определенной границе байт в оперативной памяти. Иными словами, выровненными считаются данные, адрес начала которых кратен определенному количеству байт. Для Intel Xeon Phi это 64 байта. Отметим, что в данном случае это и размер векторного регистра, и размер линейки кэшей L1 и L2.

Для определения выровненных статических массивов можно использовать следующий синтаксис [5.8]:

__declspec(align(64)) float A[1000];         //Windows
float A[1000] __attribute__((aligned(64)));  //Linux, Mac

Для выделения выровненной памяти динамически следует воспользоваться функциями _mm_malloc() и _mm_free() вместо обычных malloc() и free(). Второй аргумент функции _mm_malloc() определяет размер выравнивания в байтах:

buf = (char*) _mm_malloc(bufsize, 64);

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

Например, при передаче данных в функцию в качестве аргументов нужно сообщить компилятору о том, выровнены ли данные. Допустим, используется цикл, который обращается к массиву A как A[i], и к массиву B как B[i+n1] . Здесь i – это счетчик цикла.

Для того чтобы компилятор использовал команды работы с выровненными данными, ему необходимо сообщить, что:

Рассмотрим пример использования описанных выше конструкций:

__declspec(align(64)) float X[1000], X2[1000];

void foo(float * restrict a, int n, int n1, int n2) {
    int i;
    __assume_aligned(a, 64);
    __assume(n1%16==0);
    __assume(n2%16==0);

    for(i=0;i<n;i++) {
        X[i] += a[i] + a[i+n1] + a[i-n1]+ a[i+n2]
            + a[i-n2];
    }

    for(i=0;i<n;i++) {
        X2[i] += X[i]*a[i];
    }
}

Здесь массивы X и X2 выделены статически, поэтому дополнительно сообщать об их выравнивании не требуется. Массив a передается в функцию по указателю, а значит, компилятор ничего не знает о том, выровнен он или нет. О том, что он выровнен по границе в 64 байта мы и говорим компилятору. Дополнительно сообщаем и о выравнивании адресов со смещениями n1 и n2.

Обратите внимание на ключевое слово restrict в объявлении параметра a. Оно говорит компилятору о том, что обращение к данным по адресу a идет без каких либо перекрытий с данными из других массивов. Это позволяет компилятору определить отсутствие зависимостей в циклах.

Оба представленных в примере цикла будут автоматически векторизованы с использованием "выровненных" операций чтения/записи данных.

В качестве альтернативы предложенным конструкциям можно использовать директиву #pragma vector align перед телом векторизуемого цикла:

void test(float * a, float * b, float * c, int n)
{
  #pragma simd
  #pragma vector aligned
  for (int i = 0; i < n; i++)
    c[i] = a[i] * b[i] + a[i];
}

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

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

Пример решения этой задачи приведен ниже. Вариант кода без векторизации следующий:

static double * a;
a = _mm_malloc((N+OFFSET) * sizeof(double),64); 

#pragma omp parallel for
for (j = 0; j < N; j++)
    a[j] = 2.0E0 * a[j];

С векторизацией код будет таким:

static double * a;
a = _mm_malloc((N+OFFSET) * sizeof(double),64);

#pragma omp parallel
{
    #pragma omp master
    {
        num_threads = omp_get_num_threads();
        N1 = ((N / num_threads)/8) * num_threads * 8;
    }
}

#pragma omp parallel for
#pragma vector aligned
for (j = 0; j < N1; j++)
    a[j] = 2.0E0 * a[j];

for (j = N1; j < N; j++)
    a[j] = 2.0E0 * a[j];

Сначала мы определяем используемое количество потоков, на основании которого вычисляем новый размер цикла, позволяющий выделить каждому потоку участок с нужным выравниванием и размером.

Затем распараллеливаем и векторизуем цикл меньшего размера стандартным образом (в данном случае векторизация будет выполнена автоматически). И в завершение досчитываем оставшиеся итерации.

Подробнее о выравнивании данных для эффективной векторизации можно узнать по ссылке [5.8].

Векторизация внешних циклов

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

Векторизация внешнего цикла может быть выполнена одним из перечисленных ниже способов:

Рассмотрим следующий пример [5.10]:

for(LocalOrdinalType row=ib*BLOCKSIZE; row<top; ++row) {
    local_y[row]=0.0;

    for(LocalOrdinalType i=Arowoffsets[row];
            i<Arowoffsets[row+1]; ++i) {
        local_y[row] += Acoefs[i]*local_x[Acols[i]];
    }
}

Выполним векторизацию данного кода вторым из перечисленных выше способов:

#pragma simd
for(LocalOrdinalType row=ib*BLOCKSIZE; row<top; ++row) {
    local_y[row]=0.0;
    Inner_loop_elem_function(local_y, row, Acoefs, local_x,
        Acols, Arowoffsets);
}

__declspec(vector(uniform(Arowoffsets, Acoefs, local_x,
        Acols, local_y), linear(row))))
    Inner_loop_elem_function(float *local_y, int row, 
        float *Acoefs, float *local_x, int *Acols,
        int *Arowoffsets)
{
    for(LocalOrdinalType i=Arowoffsets[row];
            i<Arowoffsets[row+1]; ++i) {
        local_y[row] += Acoefs[i]*local_x[Acols[i]];
    }
}

Идея векторизации здесь следующая. Во-первых, описываем внутренний цикл как тело элементарной функции. Во-вторых, помещаем вызов полученной функции в тело внешнего цикла. И в-третьих, применяем simd директиву для векторизации внешнего цикла.

Подробнее о способах векторизации внешних циклов можно узнать по ссылкам [5.9,5.10].

Подходы к оптимизации прикладных программ для Intel Xeon Phi

Сопроцессоры архитектуры Intel Many Integrated Core являются мощным средством для ускорения существующих приложений, позволяя достаточно просто выполнить первичный перенос кода на сопроцессор. Однако для получения существенного прироста производительности требуется выполнить ряд шагов по оптимизации вашей параллельной программы, причем эти шаги одинаково применимы при работе как с сопроцессором, так и с обычным многоядерным процессором.

По большей части оптимизация для Intel Xeon Phi заключается в обеспечении эффективной работы с памятью.

Процесс оптимизации следует начать с выявления "узких мест" в вашей программе, для чего лучше всего воспользоваться инструментом Intel ® VTune ™ Amplifier XE.

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

Заметим, что встроенный в компилятор профилировщик циклов может использоваться только для последовательного кода. Для параллельной программы следует применять более мощные инструменты (например, Intel ® VTune ™ Amplifier XE), или отключать параллелизм.

Оптимизация приложений – это итеративный процесс, основные шаги которого представлены на рис 5.8.

Еще один подход к оптимизации заключается в поиске и устранении "узких мест" сверху вниз: от уровня операционной системы через приложение к уровню процессора (рис 5.9).

Использование встроенного профилировщика циклов

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

Для использования этой возможности необходимо включить инструментацию кода в процессе компиляции:

icc -O1 -profile-functions -profile-loops=all \
-profile-loops-report=2 ...

При этом компилятор добавит код для сбора статистики в начало и конец циклов и функций.

Далее нужно запустить приложение на интересующем вас тесте. Программа соберет статистику о своей работе и запишет ее в формате таблицы (понятном для пользователя) и в формате xml (который можно проанализировать с использованием специального GUI приложения – Loop Profile Viewer).

Файлы содержат такую информацию, как:

Пример работы GUI приложения для просмотра статистики приведен на рис 5.10.

Использование отчетов компилятора

Еще одним инструментом, полезным при оптимизации кода для процессоров и сопроцессоров Intel, являются различные отчеты компилятора.

Рассмотрим наиболее полезные в этом плане опции компиляции:

icc example.c –mmic –O3 –vec-report2

Балансировка нагрузки

Одним из существенных моментов, на которые следует обратить внимание при оптимизации приложений для Intel Xeon Phi, является балансировка нагрузки. В рамках данного раздела приведены только некоторые рекомендации касательно эффективного распределения логических потоков по вычислительным модулям сопроцессора.

Сопроцессор Intel Xeon Phi позволяет запускать одновременно 4 логических потока на ядро. Однако часто бывает эффективнее запускать меньшее их количество, т.к.:

Существуют доводы и за использование 4-х потоков на ядро:

Таким образом:

Отметим, что в режиме offload OpenMP обычно не использует ядро с номером 0 в целях повышения производительности, т.к. на нем работает операционная система и различные сервисы.

В случае, когда количество потоков вашего приложения меньше количества возможных потоков сопроцессора, следует задуматься об их распределении по ядрам. Задать желаемый алгоритм распределения средствами OpenMP можно, используя переменную окружения KMP_AFFINITY.

Возможные варианты:

Дополнительные рекомендации

Дополнительные рекомендации по оптимизации приложений для Intel Xeon Phi касаются эффективной работы с памятью:

Дополнительную информацию о программировании и оптимизации приложений для Intel Xeon Phi можно найти по ссылке [5.16].