Интернет Университет Суперкомпьютерных технологий Параллельные алгоритмы Лекция 6 Параллельные алгоритмы численного интегрирования Учебный курс Введение в параллельные алгоритмы Якобовский Михаил Владимирович проф., д.ф.-м.н. Институт прикладной математики им. М.В.Келдыша РАН, Москва

Вычислить с точностью значение определенного интеграла Пусть на отрезке [A,B] задана равномерная сетка, содержащая n+1 узел: Тогда, согласно методу трапеций, можно численно найти определенный интеграл от функции на отрезке [A,B]: Будем полагать значение J найденным с точностью, если выполнено условие: Постановка задачи Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. x1x1 x n1 =Bx1x1 x 0 =A x n2 =B x 0 =A 2 из 40 Москва, 2011 г.

IntTrap01(A,B) { n=1 J 2n =(f(A)+f(B))(B-A)/2 do { J n = J 2n n=2n s=f(A)+f(B) for(i=1;i

Пример функции Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. ABNpointseps realtime, c E E E E E E E Результаты вычисления интеграла на разных отрезках [A,B] 4 из 40 Москва, 2011 г.

main() { J= IntTrap03( A, B, f(A),f(B) ) } Адаптивный алгоритм Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. IntTrap03(A,B,fA,fB) { J=0 C=(A+B)/2 fC=f(C) sAB=(fA+fB)*(B-A)/2 sAC=(fA+fC)*(C-A)/2 sCB=(fC+fB)*(B-C)/2 sACB=sAC+sCB if(|sAB-sACB| |sACB|) J=IntTrap03(A,C,fA,fC)+IntTrap03(C,B,fC,fB) else J=sACB return J } x 2 =B x 0 =A x 1 =C Недостаток: координаты концов отрезков хранятся в программном стеке процесса и фактически недоступны программисту Преимущества: - нет повторных вычислений функции - малое число вычислений на гладких участках 5 из 40 Москва, 2011 г.

IntTrap04(A,B) { J=0 fA=f(A) fB=f(B) sAB=(fA+fB)*(B-A)/2 while(1) { Тело цикла } return J } Метод локального стека Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. C=(A+B)/2 fC=f(C) sAC=(fA+fC)*(C-A)/2 sCB=(fC+fB)*(B-C)/2 sACB=sAC+sCB if(|sAB-sACB| |sACB|) { PUT_INTO_STACK( A, C, fA, fC, sAC) A=C fA=fC sAB=sCB } else { J+=sACB if(STACK_IS_FREE) break GET_FROM_STACK( A, B, fA, fB, sAB) } B A C BA 6 из 40 Москва, 2011 г.

// данные, описывающие стек // указатель вершины стека sp=0 // массив структур в которых // хранятся отложенные задания struct { A,B,fA,fB,s } stk[1000] Процедуры и данные обеспечивающие стек Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. // макроопределения доступа к стеку #define STACK_IS_FREE (sp==0) #define PUT_INTO_STACK(A,B,fA,fB,s) { stk[sp].A=A stk[sp].B=B stk[sp].fA=fA stk[sp].fB=fB stk[sp].s=s sp++ } #define GET_FROM_STACK(A,B,fA,fB,s) { sp-- A=stk[sp].A B=stk[sp].B fA=stk[sp].fA fB=stk[sp].fB s=stk[sp].s } 7 из 40 Москва, 2011 г.

Тестирование показало, что при расчете с помощью алгоритма локального стека IntTrap04 время работы было меньше, примерно на 5%, чем при использовании IntTrap03 Примем алгоритм IntTrap04 за «наилучший» последовательный К вопросу о времени выполнения Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 8 из 40 Москва, 2011 г.

Метод геометрического параллелизма? Метод коллективного решения? ? Параллельный алгоритм интегрирования Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 9 из 40 Москва, 2011 г.

Метод геометрического параллелизма Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. main() { … for(i=0;i

p, (B-C)/ (C-A) интервал 1интервал2время1, свремя2, с 10[1e-5, ] [ , 1] [1e-5, ][ , 1] [1e-5, ][ , 1] [1e-5, ] [ , 1] [1e-5, ] [ , 1] Расчет интеграла на разных отрезках Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. B A=1e-5 C Недостаток: Большой дисбаланс загрузки процессоров 11 из 40 Москва, 2011 г.

main() { // Порождение p параллельных процессов, // каждый из которых выполняет процедуру slave for(k=0;k

Практически непригодны для решения поставленной задачи методы геометрического параллелизма (статическая балансировка) и коллективного решения (динамическая балансировка) Вывод Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 13 из 40 Москва, 2011 г.

Вычислительные системы с общей памятью Динамическая балансировка загрузки Отсутствие централизованного управления Метод глобального стека Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 14 из 40 Москва, 2011 г.

Глобальный стек Локальные стеки Стеки алгоритма Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 15 из 40 Москва, 2011 г. IntTrap04p IntTrap04p

Глобальный стек Локальные стеки Стеки алгоритма Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. Только стек Никакого процесса нет 16 из 40 Москва, 2011 г. IntTrap04p IntTrap04p

Глобальный стек Локальные стеки Стеки алгоритма Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 17 из 40 Москва, 2011 г. IntTrap04p IntTrap04p

Пусть есть доступный всем параллельным процессам список отрезков интегрирования, организованный в виде стека. Назовем его глобальным стеком. Пусть у каждого процесса есть свой, доступный только этому процессу локальный стек Перед запуском параллельных процессов в глобальный стек помещается единственная запись (в дальнейшем "отрезок"): –координаты концов отрезка интегрирования, –значения функции на концах, –приближенное значение интеграла на этом отрезке. Тогда, предлагается следующая схема алгоритма, выполняемого каждым из параллельных процессов: Пока в глобальном стеке есть отрезки: - взять один отрезок из глобального стека - выполнить алгоритм локального стека, но, если при записи в локальный стек в нем есть несколько отрезков, а в глобальном стеке отрезки отсутствуют, то: - переместить часть отрезков из локального стека в глобальный стек. Идея алгоритма Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 18 из 40 Москва, 2011 г.

- какую часть отрезков следует перемещать из локального стека в глобальный стек? - в какой момент интеграл вычислен? - что должен делать процесс у которого пуст локальный стек, если глобальный стек тоже пуст? -должен ли процесс закончить работу, если в его локальном и в глобальном стеке отрезков нет? Вопросы Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 19 из 40 Москва, 2011 г.

IntTrap04p() { // цикла обработки стека отрезков while(sdat.ntask>0) { // чтение одного интервала из списка интервалов sdat.ntask-- // указатель глобального стека GET_FROM_GLOBAL_STACK[sdat.ntask](a,b,fa,fb,sab) ИНТЕГРИРОВАНИЕ ОДНОГО ОТРЕЗКА } sdat.s_all = sdat.s_all + s } Схема Интегрирующего процесса Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 20 из 40 Москва, 2011 г.

main() {. Sem_init(sdat.sem_sum,1) //доступ к глобальной сумме открыт. } IntTrap04p() { … // Начало критической секции сложения частичных сумм sem_wait(sdat.sem_sum) sdat.s_all = sdat.s_all + s sem_post(sdat.sem_sum) // Конец критической секции сложения частичных сумм } Правильное определение общей суммы Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 21 из 40 Москва, 2011 г.

IntTrap04p() { // цикла обработки стека отрезков while(sdat.ntask>0) { // чтение одного интервала из списка интервалов sdat.ntask-- // указатель глобального стека GET_FROM_GLOBAL_STACK[sdat.ntask](a,b,fa,fb,sab) ИНТЕГРИРОВАНИЕ ОДНОГО ОТРЕЗКА } sem_wait(sdat.sem_sum) sdat.s_all = sdat.s_all + s sem_post(sdat.sem_sum) } Схема Интегрирующего процесса Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 22 из 40 Москва, 2011 г.

Схема интегрирования отрезка Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. while(1) // интегрирование одного отрезка { c=(a+b)/2;fc=f(c) sac=(fa+fc)*(c-a)/2 scb=(fc+fb)*(b-c)/2 sacb=sac+scb if(!BreakCond(sacb,sab)) // Точность на части отрезка достигнута { s+=sacb if(sp==0) break; // локальный стек пуст, выход sp--; GET_FROM_LOCAL_STACK[sp]( a, b, fa, fb, sab) } else { PUT_INTO_LOCAL_STACK[sp]( a, c, fa, fc, sac); sp++ a=c fa=fc sab=scb } if((sp>SPK) && (!sdat.ntask)) // перемещение части локального стека в общий список интервалов { while((sp>1) && (sdat.ntask

Схема интегрирования отрезка Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. c=(a+b)/2; fc=f(c) sac=(fa+fc)*(c-a)/2 scb=(fc+fb)*(b-c)/2 sacb=sac+scb // интегрирование одного отрезка while(1) { Инициализация Точность на части отрезка достигнута? Добавлять отрезки в глобальный стек? } 24 из 40 Москва, 2011 г.

Схема интегрирования отрезка Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. if(!BreakCond(sacb,sab)) { // Точность на части отрезка достигнута s+=sacb if(sp==0) break; // локальный стек пуст, выход sp--; GET_FROM_LOCAL_STACK[sp]( a, b, fa, fb, sab) } else { PUT_INTO_LOCAL_STACK[sp]( a, c, fa, fc, sac); sp++ a=c fa=fc sab=scb } // интегрирование одного отрезка while(1) { Инициализация Точность на части отрезка достигнута? Добавлять отрезки в глобальный стек? } 25 из 40 Москва, 2011 г.

Схема интегрирования отрезка Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. if((sp>SPK) && (!sdat.ntask)) // перемещение части локального стека в общий список интервалов { while((sp>1) && (sdat.ntask

while(1) { // Начало критической секции чтения из глобального // стека очередного интервала интегрирования sem_wait(sdat.sem_list) if(sdat.ntask 0) { sem_post(sdat.sem_list) // разрешить другим процессам // доступ к глобальному стеку break } sdat.ntask-- // указатель глобального стека GET_FROM_GLOBAL_STACK[sdat.ntask](a,b,fa,fb,sab) sem_post(sdat.sem_list) // Конец критической секции чтения из глобального // стека очередного интервала интегрирования … } Преждевременное окончание работы процесса Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 27 из 40 Москва, 2011 г.

Плохое условие выхода из цикла обработки стека интервалов Интегрирующие процессы не должны заканчивать работу до тех пор, пока все отрезки интервала интегрирования не будут полностью обработаны Преждевременное завершение работы приведет к получению верного ответа, но за большее время Преждевременный выход Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 28 из 40 Москва, 2011 г.

Отрезок интегрирования может находиться в нескольких состояниях: - находится в глобальном стеке интервалов; - обрабатывается некоторым интегрирующим процессом; - находится в локальном стеке интервалов некоторого процесса; - полностью обработан: известно значение интеграла на этом отрезке и оно прибавлено к локальной частичной сумме соответствующего процесса. "Время жизни" отрезка, после того, как некоторый процесс начал его обработку, относительно невелико - отрезок разбивается на две части и перестает существовать, породив два новых отрезка. Таким образом, требование "все отрезки интервала интегрирования полностью обработаны" означает, что: - функция проинтегрирована на всех отрезках, покрывающих исходный интервал интегрирования; - полученные на отрезках интегрирования значения интегралов добавлены к частичным суммам соответствующих процессов. Если глобальный и локальный стеки пусты Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 29 из 40 Москва, 2011 г.

семафоры доступа: sdat.sem_listсемафор доступа к глобальному стеку отрезков sdat.s_allсемафор доступа к значению интеграла семафоры состояния: sdat.sem_task_presentсемафор наличия записей в глобальном стеке отрезков переменные: sdat.list_of_tasksглобальный стек отрезков sdat.ntaskчисло записей в глобальном стеке отрезков - указатель глобального стека отрезков sdat.nactiveчисло активных процессов sdat.s_allзначение интеграла Необходимые данные Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 30 из 40 Москва, 2011 г.

int slave_thr() { // все переменные, начинающиеся с sdat. являются глобальными, // к ним имеет доступ каждый из запущенных процессов slave_thr // и запускающая программа main // sp, s - локальные переменные процессов slave_thr sp=0 // указатель локального стека - локальный стек пуст s=0 // частичная сумма интегралов, вычисленных на отрезках, // обработанных данной копией процесса // начало цикла обработки стека интервалов while(1) { // ожидание появления в глобальном стеке интервалов для обработки sem_wait(sdat.sem_task_present) // чтение одного интервала из списка интервалов // Начало критической секции чтения из глобального // стека очередного интервала интегрирования // sem_wait(sdat.sem_list) sdat.ntask-- // указатель глобального стека GET_OF_GLOBAL_STACK[sdat.ntask](a,b,fa,fb,sab) if(sdat.ntask) sem_post(&sdat.sem_task_present) if(ab) // отрезок является терминальным break // выйти из цикла обработки стека интервалов // начало цикла интегрирования одного интервала while(1) { c=(a+b)/2 fc=fun(c) sac=(fa+fc)*(c-a)/2 scb=(fc+fb)*(b-c)/2 sacb=sac+scb if(!BreakCond(sacb,sab)) { s+=sacb if(!sp) // локальный стек пуст break // выход из цикла интегрирования // одного интервала sp-- GET_FROM_LOCAL_STACK[sp]( a, b, fa, fb, sab) } else { PUT_TO_LOCAL_STACK[sp]( a, c, fa, fc, sac) sp++ a=c fa=fc sab=scb } // перемещение части локального стека // в общий список интервалов if((sp>SPK) && (!sdat.ntask)) { // Начало критической секции заполнения глобального // стека отрезками интегрирования // sem_wait(sdat.sem_list) if(!sdat.ntask) { // установить семафор наличия // записей в глобальном стеке sem_post(sdat.sem_task_present) } while((sp>1) && (sdat.ntaskb) // sem_wait(&sdat.sem_list) sdat.nactive-- if( (!sdat.nactive) && (!sdat.ntask) ) { // запись в глобальный стек списка терминальных отрезков for(i=0;i

Инициализация Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. int slave_thr() { // все переменные, начинающиеся с sdat. являются глобальными, // к ним имеет доступ каждый из запущенных процессов slave_thr // и запускающая программа main // sp, s - локальные переменные процессов slave_thr sp=0 // указатель локального стека - локальный стек пуст s=0 // частичная сумма интегралов, вычисленных на отрезках, // обработанных данной копией процесса 32 из 40 Москва, 2011 г.

Доступ процесса IntTrap04p к глобальному стеку Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. // начало цикла обработки глобального стека while(1) { // ожидание появления в глобальном стеке интервалов для обработки sem_wait(sdat.sem_task_present) // чтение одного интервала из списка интервалов sem_wait(sdat.sem_list) sdat.ntask-- // указатель глобального стека GET_OF_GLOBAL_STACK[sdat.ntask](a,b,fa,fb,sab) if(sdat.ntask) sem_post(&sdat.sem_task_present) if(ab) // отрезок является терминальным break // выйти из цикла обработки стека интервалов 33 из 40 Москва, 2011 г.

Запись терминальных отрезков Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. // Начало критической секции заполнения глобального // стека терминальными отрезками (a>b) sem_wait(&sdat.sem_list) sdat.nactive-- if( (!sdat.nactive) && (!sdat.ntask) ) { // запись в глобальный стек списка терминальных отрезков for(i=0;i

Время выполнения Ускорение Эффективность Np1234 tiger.jscc.ru ga03.imamod.ru Результаты тестирования Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. Np1234 tiger.jscc.ru ga03.imamod.ru np1234 tiger.jscc.ru 100%101%102%100% ga03.imamod.ru 100%99% 35 из 40 Москва, 2011 г.

- какую часть отрезков следует перемещать из локального стека в глобальный стек? - в какой момент интеграл вычислен? - что должен делать процесс у которого пуст локальный стек, если глобальный стек тоже пуст? -должен ли процесс закончить работу, если в его локальном и в глобальном стеке отрезков нет? - На сколько отрезков изначально разбить отрезок? Вопросы Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 36 из 40 Москва, 2011 г.

Рассмотрен ряд методов вычисления интегралов на многопроцессорных системах, проанализированы их преимущества и недостатки Показано, что методы геометрического параллелизма и коллективного решения неприменимы для эффективного численного интегрирования функций общего вида Показано, что идея децентрализованного управления расчетом эффективно решает задачу Заключение Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 37 из 40 Москва, 2011 г.

Якобовский М.В., Кулькова Е.Ю. Решение задач на многопроцессорных вычислительных системах с разделяемой памятью. - М.: СТАНКИН, – 30 c. Литература Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 38 из 40 Москва, 2011 г.

Якобовский М.В. проф., д.ф.-м.н., зав. сектором «Программного обеспечения многопроцессорных систем и вычислительных сетей» Института прикладной математики им. М.В.Келдыша Российской академии наук mail: web: Контакты Введение в параллельные алгоритмы: Параллельные алгоритмы интегрирования функций © Якобовский М.В. 39 из 40 Москва, 2011 г.