Лекция 2: Принципы переноса прикладных программных пакетов на Intel Xeon Phi

Прямой перенос программного пакета на сопроцессор

Кратко опишем особенности исходной версии программного пакета.

Распараллеливание ведется по фотонам, каждый поток обсчитывает свой набор траекторий. Фотоны делятся поровну между потоками:

void LaunchOMP(InputInfo* input, OutputInfo* output,
    MCG59* randomGenerator, int numThreads)
{
    omp_set_num_threads(numThreads);

    …

    #pragma omp parallel
    {
        int threadId = omp_get_thread_num();

        for (uint64 i = threadId;
            i < input->numberOfPhotons; i += numThreads)
        {
            ComputePhoton(specularReflectance, input,
                &(threadOutputs[threadId]),
                &(randomGenerator[threadId]), trajectory);
        }
    }
    …
}

Работа ведется с числами двойной точности.

Результатами трассировки являются:

• Величина сигнала на детекторе – массив, количество элементов которого равно числу детекторов, размер для 10 детекторов – 80 Б.
• Общая фотонная карта траекторий – массив, количество элементов которого равно числу ячеек сетки. Для хранения результатов используется трехмерная сетка, сохраняются данные о количестве посещений фотонами каждой из ее ячеек ( рис. 2.10). Если рассматривать сетку размером 100*100*50 элементов, размер ее будет равен 4 МБ.
• Фотонные карты траекторий для каждого детектора. Фотонная карта для каждого детектора хранится в массиве, аналогичном массиву с общей фотонной картой. Размер результатов для 10 детекторов – 40 МБ.

увеличить изображение
Рис. 2.10. Хранение числа посещений фотоном ячеек сетки в самой сетке

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

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

И если для исполнения на центральном процессоре объем дополнительной памяти будет равен 350 МБ при использовании 8 потоков, то для работы на сопроцессоре Intel Xeon Phi с 240 потоками понадобится хранить уже более 10 ГБ данных.

Еще одна особенность алгоритма состоит в необходимости хранения траектории каждого отдельного фотона в процессе его трассировки. Это связано с тем, что до окончания трассировки фотона узнать, попал ли он в детектор, нельзя. Соответственно, если фотон попал в детектор, траектория его движения суммируется с общей картой траекторий для данного детектора, а если нет – игнорируется.

В исходной версии выделяется массив для хранения траектории для каждого потока. С точки зрения структур данных используется все та же трехмерная сетка, соответственно размер массива равен 4 МБ. Перед началом трассировки каждого отдельного фотона этот массив обнуляется.

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

Тестовые запуски проводились на следующем оборудовании: CPU Intel Xeon E5-2690 2.9 ГГц (8 ядер), Intel Xeon Phi 7110X, 64 ГБ RAM, использовался компилятор Intel C/C++ Compiler 14.0, операционная система – CentOS 6.2.

Во всех экспериментах расчеты выполнялись с двойной точностью.

Результаты прямого переноса приведены на рис. 2.11.

увеличить изображение
Рис. 2.11. Время работы исходной версии программы на CPU и Intel Xeon Phi

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

Оптимизация 1: обновление структуры данных для хранения траектории фотона в процессе трассировки

Перед началом оптимизации имеет смысл выявить наиболее медленные участки программы. Для этого воспользуемся профилировщиком Intel VTune Amplifier XE. Результаты профилировки приведены на рис. 2.12.

Рис. 2.12. Результаты профилировки базовой версии программы

Как видно из графиков, наибольшее время занимает функция memset(), которая используется в момент начала процесса трассировки фотона для обнуления памяти, предназначенной для хранения траектории его движения. Напомним, что текущая траектория движения фотона хранится в виде трехмерной сетки ( рис. 2.10).

увеличить изображение
Рис. 2.13. Хранение числа посещений фотоном ячеек сетки в виде списка координат этих ячеек

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

Для устранения описанных выше недостатков предлагается использовать другую структуру данных для хранения траектории фотона: следует хранить не число посещений данной ячейки для всей сетки, а список координат посещенных ячеек ( рис. 2.13).

Экспериментальные данные показывают, что максимальная длина траектории фотона в данном примере не превосходит 2048 шагов (эта величина существенно зависит от параметров биотканей и размера сетки). Размер структуры данных в этом случае составляет всего 6 КБ.

Результаты применения этой оптимизации приведены на рис. 2.14.

увеличить изображение
Рис. 2.14. Время работы программы после оптимизации структуры данных для хранения траектории движения фотона на CPU и Intel Xeon Phi

Изменение используемой структуры данных позволило сократить время вычислений на CPU в 1 поток на 18%, а в 8 – на 9%. Время работы программы на MIC в 120 потоков сократилось на 31%.