Параллельные методы решения дифференциальных уравнений в частных производных

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

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

T1=kmN2,

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

11.2.1. Использование OpenMP для организации параллелизма

Рассмотрим возможные способы организации параллельных вычислений для сеточных методов на многопроцессорных вычислительных системах с общей памятью. При изложении материала будем предполагать, что имеющиеся в составе системы процессоры обладают равной производительностью, являются равноправными при доступе к общей памяти и время доступа к памяти является одинаковым (при одновременном доступе нескольких процессоров к одному и тому же элементу памяти очередность и синхронизация доступа обеспечиваются на аппаратном уровне). Как уже отмечалось ранее, многопроцессорные системы подобного типа обычно именуются симметричными мультипроцессорами ( symmetric multiprocessors, SMP ) – см. п. 1.3.1.

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

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

Рассмотренный выше подход является основой технологии OpenMP (см., например, [ [ 27 ] ]), наиболее широко применяемой в настоящее время для организации параллельных вычислений на многопроцессорных системах с общей памятью. В рамках данной технологии директивы параллелизма используются для выделения в программе параллельных областей ( parallel regions ), в которых последовательный исполняемый код может быть разделен на несколько раздельных командных потоков ( threads ). Далее эти потоки могут исполняться на разных процессорах вычислительной системы. В результате такого подхода программа представляется в виде набора последовательных ( однопотоковых ) и параллельных ( многопотоковых ) участков программного кода (см. рис. 11.3). Подобный принцип организации параллелизма получил наименование "вилочного" ( fork-join ) или пульсирующего параллелизма. Более полная информация по технологии OpenMP может быть получена в литературе (см., например, [ [ 27 ] , [ 66 ] ]) или в информационных ресурсах сети Интернет. В данной лекции возможности OpenMP будут излагаться в объеме, необходимом для демонстрации возможных способов разработки параллельных программ для рассматриваемого учебного примера решения задачи Дирихле.

11.2.2. Проблема синхронизации параллельных вычислений

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

Алгоритм 11.2. Первый вариант параллельного алгоритма Гаусса – Зейделя

// Алгоритм 11.2
omp_lock_t dmax_lock;
omp_init_lock (dmax_lock);
do {	
  dmax = 0; // максимальное изменение значений u
#pragma omp parallel for shared(u,N,dmax) private(i,temp,d)
  for ( i=1; i<N+1; i++ ) {
#pragma omp parallel for shared(u,N,dmax) private(j,temp,d)
    for ( j=1; j<N+1; j++ ) {
      temp = u[i][j];
      u[i][j] = 0.25*(u[i-1][j]+u[i+1][j]+
                u[i][j-1]+u[i][j+1]–h*h*f[i][j]);
      d = fabs(temp-u[i][j])
      omp_set_lock(dmax_lock);
        if ( dmax < d ) dmax = d;
      omp_unset_lock(dmax_lock);
    } // конец вложенной параллельной области
  } // конец внешней параллельной области
} while ( dmax > eps );

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

Рис. 11.3. Параллельные области, создаваемые директивами OpenMP

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

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

Результаты вычислительных экспериментов приведены в табл. 11.1 (здесь и далее для параллельных программ, разработанных с использованием технологии OpenMP, использовался четырехпроцессорный сервер кластера Нижегородского университета с процессорами Pentium III, 700 Mhz, 512 RAM).

Оценим полученный результат. Разработанный параллельный алгоритм является корректным, т.е. обеспечивающим решение поставленной задачи. Использованный при разработке подход обеспечивает достижение практически максимально возможного параллелизма – для выполнения программы может быть задействовано вплоть до N² процессоров. Тем не менее результат не может быть признан удовлетворительным: программа будет работать медленно и вместо ускорения мы получим замедление вычислений. Основная причина такого положения дел – чрезмерно высокая синхронизация параллельных участков программы. В нашем примере каждый параллельный поток после усреднения значений u_ij должен проверить (и возможно, изменить) значение величины dmax. Разрешение на использование переменной может получить только один поток – все остальные потоки должны быть блокированы. После освобождения общей переменной управление может получить следующий поток и т.д. В результате необходимости синхронизации доступа многопотоковая параллельная программа превращается фактически в последовательно выполняемый код, причем менее эффективный, чем исходный последовательный вариант, т.к. организация синхронизации приводит к дополнительным вычислительным затратам – см. рис. 11.4. Следует обратить внимание, что, несмотря на идеальное распределение вычислительной нагрузки между процессорами, для приведенного на рис. 11.4 соотношения параллельных и последовательных вычислений, в каждый текущий момент времени (после момента первой синхронизации) только не более двух процессоров одновременно выполняют действия, связанные с решением задачи. Подобный эффект вырождения параллелизма из-за интенсивной синхронизации параллельных участков программы обычно именуется сериализацией ( serialization ).

( k -количество итераций, t -время (сек), S - ускорение )

Рис. 11.4. Пример возможной схемы выполнения параллельных потоков при наличии синхронизации (взаимоисключения)

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

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

Алгоритм 11.3. Второй вариант параллельного алгоритма Гаусса – Зейделя

// Алгоритм 11.3
omp_lock_t dmax_lock;
omp_init_lock(dmax_lock);
do {
  dmax = 0; // максимальное изменение значений u
#pragma omp parallel for shared(u,N,dmax) private(i,temp,d,dm)
  for ( i=1; i<N+1; i++ ) {
    dm = 0;
    for ( j=1; j<N+1; j++ ) {
      temp = u[i][j];
      u[i][j] = 0.25*(u[i-1][j]+u[i+1][j]+
                u[i][j-1]+u[i][j+1]–h*h*f[i][j]);
      d = fabs(temp-u[i][j])
      if ( dm < d ) dm = d;
    }
    omp_set_lock(dmax_lock);
      if ( dmax < dm ) dmax = dm;
    omp_unset_lock(dmax_lock);
    } 
  } // конец параллельной области
} while ( dmax > eps );

Как результат выполненного изменения схемы вычислений, количество обращений к общей переменной dmax уменьшается с N² до N раз, что должно приводить к существенному снижению затрат на синхронизацию потоков и уменьшению проявления эффекта сериализации вычислений. Результаты экспериментов с данным вариантом параллельного алгоритма, приведенные в табл. 11.1, показывают существенное изменение ситуации – ускорение в ряде экспериментов оказывается даже большим, чем используемое количество процессоров (такой эффект сверхлинейного ускорения достигается за счет наличия у каждого из процессоров вычислительного сервера своей быстрой кэш-памяти). Следует также обратить внимание, что улучшение показателей параллельного алгоритма достигнуто при снижении максимально возможного параллелизма (для выполнения программы может использоваться не более N процессоров).