Участник:Bannikov Nikolay/Параллельный алгоритм вычисления всех частных сумм

Свойства и структура алгоритма

Общее описание алгоритма

Данный алгоритм предназначен для вычисления всех частных сумм заданного массива данных. Для входного массива:

[math]\displaystyle{ X = \{x_1,x_2,...,x_n\} }[/math]

требуется найти массив:

[math]\displaystyle{ S = \sum_{i=1}^n x_i }[/math]

такой, что каждый его элемент [math]\displaystyle{ S[i] }[/math] равен сумме всех элементов [math]\displaystyle{ X }[/math] от первого до [math]\displaystyle{ i }[/math]-го включительно.
Алгоритм относится к классу алгоритмов сдваивания. Его ключевая идея заключается в том, что на каждом шаге расстояние, на которое передаются и суммируются уже накопленные частичные результаты, удваиваются.
Входные данные представляют собой массив числовых элементов, для которых определена операция сложения, однако алгоритм может быть обобщен на любую ассоциативную бинарную операцию (умножение, нахождение минимума/максимума, логическое И/ИЛИ и.т.д.).
Алгоритм широко используется в параллельных вычислениях.

Математическое описание алгоритма

Пусть дан массив чисел:

[math]\displaystyle{ X = \{x_1,x_2,...,x_n\} }[/math]

Необходимо найти массив

[math]\displaystyle{ S = \{s_1,s_2,...,s_n\} }[/math]

где [math]\displaystyle{ s_i = x_1 + x_2 + ... + x_i }[/math] для всех [math]\displaystyle{ i }[/math] от 1 до n.

Алгоритм выполняется за [math]\displaystyle{ k = log_{2} n }[/math] шагов:

1. Каждый процессор [math]\displaystyle{ i }[/math] хранит значение [math]\displaystyle{ s_i = x_i }[/math]
2. Для каждого шага [math]\displaystyle{ t }[/math] от 1 до k:
   • Вычисляем расстояние [math]\displaystyle{ d = 2^{(t-1)} }[/math]
   • Для каждого индекса [math]\displaystyle{ i \gt d }[/math] процессор [math]\displaystyle{ i }[/math] выполняет: [math]\displaystyle{ s_i = s_i + s_{i-d} }[/math]
   • Для каждого индекса [math]\displaystyle{ i ≤ d }[/math] значение [math]\displaystyle{ s_i }[/math] остается неизменным
3. После [math]\displaystyle{ k }[/math] шагов массив [math]\displaystyle{ S }[/math] содержит все частные суммы

Простой пример для [math]\displaystyle{ n = 4: }[/math]

• Начало: [math]\displaystyle{ [x_1,x_2,x_3,x_4] }[/math]
• Шаг 1 (d=1): [math]\displaystyle{ [x_1, x_1 + x_2, x_2 + x_3, x_3 + x_4] }[/math]
• Шаг 2 (d=2): [math]\displaystyle{ [x_1, x_1 + x_2, x_1 + x_2 + x_3, x_1 + x_2 + x_3 + x_4] }[/math]

Вычислительное ядро алгоритма

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

[math]\displaystyle{ s_i = s_i + s_{i-d} }[/math]

где:

• [math]\displaystyle{ s_i }[/math] - текущее значение в процессоре [math]\displaystyle{ i }[/math]
• [math]\displaystyle{ s_{i-d} }[/math] - значение из процессора, находящегося на расстоянии [math]\displaystyle{ d }[/math] слева
• [math]\displaystyle{ d }[/math] - расстояние между процессорами, которое удваивается на каждом шаге

Характеристики вычислительного ядра:

1. Простота: Каждая операция - одно сложение
2. Регулярность: На каждом шаге выполняется одинаковый тип операций, только меняется расстояние
3. Параллельность: Все операции на одном шаге могут выполняться одновременно разными процессорами
4. Коммуникационная интенсивность: Основная сложность - в организации быстрого обмена данными между процессорами, а не в самих вычислениях

На первом шаге выполняется [math]\displaystyle{ n - 1 }[/math] сложений, на втором [math]\displaystyle{ n - 2 }[/math] сложений и так далее. Несмотря на то, что общее количество операций [math]\displaystyle{ (n · log_{2} n) }[/math] больше чем в последовательном алгоритме [math]\displaystyle{ n }[/math] операций, время выполнения сокращается с [math]\displaystyle{ n }[/math] шагов до [math]\displaystyle{ log_{2} n }[/math] шагов благодаря параллелизму.

Макроструктура алгоритма

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

Алгоритм представляет собой последовательность [math]\displaystyle{ log_{2} n }[/math] итераций.
Каждая итерация для всех активных процессоров представима следующими двумя этапами:

1. Передача данных на расстояние [math]\displaystyle{ d }[/math]
2. Сложение: [math]\displaystyle{ s[i] = s[i] + }[/math] полученное_значение

Схема реализации последовательного алгоритма

Для заданного массива чисел [math]\displaystyle{ X = [x_1, x_2,...,x_n] }[/math] требуется вычислить массив частичных сумм [math]\displaystyle{ S = [s_1,s_2,...,s_n] }[/math], где

[math]\displaystyle{ S_k = \sum_{i=1}^n x_i, 1 ≤ k ≤ n }[/math]

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

[math]\displaystyle{ S_k = S_{k-1} + x_k, \quad k = 1,2,\ldots,n; \quad S_0 = 0 }[/math]

Основной вычислительный фрагмент кода на С++:

    // Вычисление частичных сумм
    s[0] = a[0];
    for (int i = 1; i < n; i++) 
    {
        s[i] = s[i - 1] + a[i];
    }

• Последовательный алгоритм выполняет [math]\displaystyle{ O(n) }[/math] операций за [math]\displaystyle{ O(n) }[/math] времени на 1 процессоре
• Параллельный алгоритм выполняет [math]\displaystyle{ O(n \cdot \log n) }[/math] операций за [math]\displaystyle{ O(\log n) }[/math] времени на [math]\displaystyle{ n }[/math] процессорах
• При [math]\displaystyle{ n = 1 }[/math] или при выполнении на одном ядре параллельный алгоритм вырождается в последовательный

Рис. 1 Диаграмма маршрутизации последовательного алгоритма вычисления всех частных сумм

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

Последовательная сложность алгоритма

Последовательная сложность алгоритма вычисления всех частичных сумм линейна по размеру входных данных: требуется [math]\displaystyle{ n - 1 }[/math] операции сложения, где [math]\displaystyle{ n }[/math] - размер массива и [math]\displaystyle{ O(n) }[/math] обращений к памяти.

Информационный граф

Рис. 2 Диаграмма маршрутизации параллельного алгоритма вычисления всех частных сумм

На рис.2 изображён граф алгоритма. В данном случае выполнено суммирование 8 элементов массива.

• Видно, что параллельный алгоритм вычисления всех частных сумм требует [math]\displaystyle{ log_{2} n }[/math] сложений со степенью параллелизма [math]\displaystyle{ n }[/math], если для обеспечения однородности вычислительной схемы на ПЭ (процессорном элементе), для которых нет данных слева, отправлять нулевые значения.
• Каждая операция алгоритма имеет максимально возможную степень параллелизма [math]\displaystyle{ n }[/math]
• Общее число скалярных арифметических операций увеличилось с [math]\displaystyle{ n }[/math] до [math]\displaystyle{ (n · log_{2} n) }[/math], а время выполнения задачи уменьшилось

Таким образом, увеличение степени параллелизма на каждом этапе от 1 до n происходит за счет существенного увеличения числа скалярных арифметических операций.

Ресурс параллелизма алгоритма

Для вычисления всех частных сумм массива порядка [math]\displaystyle{ n }[/math] с помощью алгоритма в его параллельном варианте требуется выполнить [math]\displaystyle{ (\log_{2} n) }[/math] ярусов с убывающим параллелизмом: от [math]\displaystyle{ n - 1 }[/math] до [math]\displaystyle{ 1 }[/math] операции суммирования на ярус. В предположении доступности неограниченного числа необходимых процессоров требуется [math]\displaystyle{ O(n \log n) }[/math] операций при [math]\displaystyle{ p = n }[/math] за [math]\displaystyle{ O(\log n) }[/math] времени. При классификации по высоте ЯПФ алгоритм имеет логарифмическую высоту [math]\displaystyle{ O(\log n) }[/math],а при классификации по ширине ЯПФ алгоритм имеет линейную ширину ярусов [math]\displaystyle{ O(n) }[/math], реализуя массовый координатный параллелизм по индексу процессорных элементов. Максимальное ускорение ограничено [math]\displaystyle{ O(\frac{n}{\log n}) }[/math] из закона Амдала. Наиболее эффективен при [math]\displaystyle{ n ≫ p ≥ \log n }[/math]. Поскольку эффективность алгоритма уменьшается при росте [math]\displaystyle{ n }[/math] и составляет [math]\displaystyle{ O(\frac{1}{\log n}) }[/math].

Входные и выходные данные алгоритма

Входные данные

• Объект: одномерный массив чисел [math]\displaystyle{ X = [x_1,x_2,...,x_n] }[/math]
• Тип: Любой числовой (целые, вещественные, комплексные)
• Размер: [math]\displaystyle{ n }[/math] элементов ([math]\displaystyle{ n ≥ 0 }[/math])
• Объем памяти: [math]\displaystyle{ n × sizeof }[/math](тип)

Выходные данные

• Объект: Массив частичных сумм [math]\displaystyle{ S = [s_1,s_2,...,s_n] }[/math]
• Формула: [math]\displaystyle{ s_i = x_1 + x_2 + ... + x_i }[/math]
• Размер: [math]\displaystyle{ n }[/math] элементов
• Объем памяти: [math]\displaystyle{ n × sizeof }[/math](тип)

Ключевые особенности

1. Вход и выход одинакового размера и типа
2. Минимум [math]\displaystyle{ 2n }[/math] элементов памяти (вход + выход)
3. Можно перезаписать входной массив с целью экономии памяти
4. Последний элемент [math]\displaystyle{ s_n }[/math] содержит сумму всех входных чисел

Ограничения

• Переполнение: Риск для целых чисел при больших [math]\displaystyle{ n }[/math]
• Точность: Накопление ошибок для чисел с плавающей точкой
• Память: Максимальный [math]\displaystyle{ n }[/math] ограничен доступной памятью

Свойства алгоритма

Соотношение последовательной и параллельной сложности в предположении доступности неограниченного числа необходимых процессоров, является [math]\displaystyle{ \frac{n}{\log_{2} n} }[/math] отношение линейной к логарифмической, ограничено законом Амдала. Алгоритм дает существенное теоретическое ускорение, почти линейное с эффективностью [math]\displaystyle{ \frac{1}{\log_{2} n} }[/math]. При этом вычислительная мощность алгоритма, как отношение числа арифметических операций к суммарному объему входных и выходных данных растет как [math]\displaystyle{ \log n }[/math]: при больших [math]\displaystyle{ n }[/math] на один элемент данных приходится все больше операций. Это положительный момент для реализации на ускорителях и сопроцессорах: накладные расходы на пересылку данных относительно стоимости вычислений снижаются. Для целых чисел алгоритм абсолютно устойчив, получаем точный результат при отсутствии переполнения. Для чисел с плавающей точкой условно устойчив, поскольку накапливается ошибка округления, где погрешность [math]\displaystyle{ O(n · ε) }[/math], где [math]\displaystyle{ ε }[/math] - машинная точность.
Алгоритм сбалансирован, поскольку лишь одна операция сложения, нет умножений/делений и.т.д. Алгоритм полностью детерминирован, поскольку число этапов фиксировано и равно [math]\displaystyle{ \log_{2} n }[/math], набор операций, выполняемых каждым ПЭ на каждом этапе заранее известен и не зависит от входных данных, отсутствуют ветвления по данным, циклы с выходом по точности, использование случайных чисел и.т.п.
В алгоритме нет «массовых рассылок» скалярных величин с большой степенью исхода одной вершины. Это снижает требования к механизмам широковещательной передачи и упрощает эффективный доступ к данным: каждый результат либо остаётся локальным, либо передаётся ограниченному числу соседей по коммуникационной схеме.
В линейной системе координат присутствуют длинные дуги, что требует аккуратного отображения на реальную топологию, но при подходящем выборе топологии эти дуги становятся локальными.
Граф алгоритма регулярно устроен и допускает компактные укладки и конвейеризацию при реализации на специализированных архитектурах.

Литература

1. А.В.Старченко, В.Н.Берцун. Методы параллельных вычислений. Изд-во Том. ун-та, 2013. - 223 с. (Серия "Суперкомпьютерное образование")