Параллельная реализация явных методов типа Рунге-Кутта Текст научной статьи по специальности «Математика»

3. Шапиро В.Я., Григорьев И.В., Жукова А.И. Влияние сдвиговых деформаций на процесс циклического уплотнения почвы // Естественные и технические науки. - 2006. - № 1(21). - С.174-180.

4. Григорьев И.В. Снижение отрицательного воздействия на почву колесных трелевочных тракторов обоснованием режимов их движения и технологического оборудования. - СПб.: ЛТА, 2006. - 236 с.

5. Средощадящие технологии разработки лесосек в условиях Северо-Западного региона Российской Федерации / И.В. Гоигорьев [и др.]. - СПб.: ЛТА, 2008. - 176 с.

6. Рудов С.Е. Пути исследования воздействия лесозаготовительных машин на почву при несплошных рубках // Леса России в XXI веке: мат-лы 1-й междунар. науч.-практ. конф. - СПб.: ЛТА, 2009. - С. 205-211.

7. Иванов В.А., Вовченко Н.Д. Амплитудно-частотные характеристики малого лесного трактора // Леса России в XXI веке: мат-лы 1-й междунар. науч.-практ. конф. - СПб.: ЛТА, 2009. - С. 205-211.

8. Бленд Д. Теория линейной вязкоупругости. - М.: Мир, 1965. - 199 с.

УДК 519.852.6 Г.В. Ващенко, Е.А. Новиков

ПАРАЛЛЕЛЬНАЯ РЕАЛИЗАЦИЯ ЯВНЫХ МЕТОДОВ ТИПА РУНГЕ-КУТТА*

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

Ключевые слова: явные методы Рунге-Кутта, параллельные алгоритмы, кластеры, полный граф.

G.V. Vashchenko, Ye.A. Novikov PАRALLEL IMPLEMENTATION OF THE EXPLICIT TECHNIQUES OF RUNGE-KUTT TYPE

Explicit techniques of Runge-Kutt type for the numerical initial problem solution for the systems of the ordinary differential equations of the first order are considered in the article. The parallel computing schemes of techniques directed to the multiprocessing computing systems of cluster architecture are offered.

Key words: explicit Runge-Kutt techniques, parallel algorithms, clusters, completely-connected graph.

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

У = Я^уХ у(*о)=у°- (1)

Для численного решения задачи (1) применяются явные в-стадийные методы типа Рунге-Кутта, (п+1)-й шаг в которых задается формулами

г = 1

г-1

*‘"’=/(',+<А+,,У"’+А,+.1л*Г), (2)

7-1

1=1, 2,..., я.

Конкретный метод Рунге-Кутта описывается набором коэффициентов Ь, с, ац, 1<^, 2</<(/-1) [2-4].

'Работа поддержана грантами РФФИ №08-01-00621 и Президента НШ-3431.2008.9.

Здесь представлены параллельные вычислительные схемы для в-стадийных явных методов типа Рунге-Кутта и параллельные алгоритмы, ориентированные на применение в многопроцессорных вычислительных системах кластерной архитектуры. Основной подход при конструировании параллельных методов состоял в распараллеливании последовательных численных алгоритмов, использовании декомпозиции на подзадачи [1; 5] исходной задачи (2) и анализе информационных взаимосвязей между ними. Поиск максимально возможного независимого набора операций с целью их распараллеливания осуществлялся путем построения графов зависимостей. Суть построения графов заключается в сопоставлении вершинам графа выполняемых независимо операций, а дугам - результата влияния выполнения операции на результат некоторых соседних операций.

Последовательная вычислительная схема и алгоритм. В данной работе будем отталкиваться от неоднократно апробированных последовательных алгоритмов решения задачи (1) методами (2). Рассмотрим последовательную вычислительную схему. Для определенности зададимся некоторым отрезком [йз, Т] и введем равномерную сетку \/\/п=(Ь, и, й„) с величиной шага Лп+1=/7. На сетке и/„ в начальный момент вре-

мени

йо в качестве начального условия зададим вектор у = (у{ , у2 , ..ук ). Тогда построение численного решения сводится к нахождению последовательности значений у(0), у(1),у(п), являющихся приближением значений у(ґ0 ), у(^ ), у(^2 )., у(їп ) точного решения у(й) в узлах сетки

Определение значений компонент вектора приближенного решения

~(я+1) /г;(и+1) С,(и+1) ~(к!+1) л о. “ “

У ~(X ■•У2 ?---?Уы ) осуществляется по формуле (2), записанной для вычисления каждой

компоненты вектора у , то есть

Я'"1’ = У;>+1г\Ь,К1"]+Ь2К?>+...+Ь1К^1

(3)

УТ=УМо), У=1,2,...,ДГ,

где

4"]=/, (1Лс2Кум+а2_ №<">);

£<">=/(? +сИ,у(’) + Ыа Ж\п)+а ,£<”>+...+а ,К(’] ).

5,7 *1 7^ п 5 V. 5,1 1,7 5,2 2,7 5,5-1 5-1,7/

Как видно из (2) и (3), определение значения у{п+1) сводится к строго последовательному вычислению коэффициентов К[п),К2и),..Кп), их умножению на коэффициенты Ь, /=1, 2, ..., 5 соответственно

и последующему суммированию.

Исходя из формулы (3), последовательный алгоритм в общем случае построения численного решения задачи (1) может быть представлен в виде:

начало;

ввод начального вектора у0; вычисление шага Л;

цикл по числу узлов сетки п = 1, 2, ..., М - 1; цикл по числу стадий т = 1, 2, ..., 5; цикл по числу компонент / = 1, 2, ..., N

вычисление коэффициентов К”) ;

конец цикла по /; конец цикла по т;

цикл по числу компонент ] = 1, 2, ..., N

, ,0+1) вычисление гЛ ; ;

конец цикла по /;

цикл по числу компонент 7 = 1,2, ..., А/;

,,(и) (и+1)

сохранение, _у\ = У) ,

конец цикла по 7;

вывод ;

конец цикла по п;

конец.

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

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

Параллельная вычислительная схема. Для разработки параллельной схемы решения задачи (1) явными методами (2) заметим, что последовательная схема (3) дает основание организации двух типов параллельных вычислений:

а) вычисление отдельных компонент коэффициентов К2<п>,..., К5(п) и вектора численного реше-

б) вычисление отдельных операций внутри одного шага метода.

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

Анализ формулы (3) показывает, что вычисление любой компоненты вектора приближенного решения

т=1, 2, ..., 5; I = 1, 2, ..., N и последующему суммированию результатов умножения ЬК1, I + ЬгК^ +...+&&,I. В силу этого, за базовую подзадачу примем вычисления компонент Кт, I, т = 1, 2, ..., 5; I = 1, 2,., N.

Введем в рассмотрение вектор

Теперь построим общую параллельную схему вычислений. Рассмотрим информационные зависимо-

ния у

,-,(«+1) _

_у(й + 1) В каждом узле сетки Жп СВОДИТСЯ К вычислению произведения коэффициента Ьт на компоненты Кп, /,

т — 1

с элементами

т-1

(4)

т = 1,2,...,5, / = 1,2,...,А/-,

тогда коэффициенты Кбудут определяться формулой

(5)

сти, возникающие при вычислении коэффициентов_К, т=1,2,...,в, /=1,2,..., N, при получении вектора у{п+1). Будем предполагать, что имеется N (р=Л/) компьютеров (рис.1), на каждом из которых вычисляется последовательно один коэффициент Кпо формуле (5) и элемент бП , т=1,2, ..., 5, / = 1, 2,., N по

формуле (4). В системе каждый компьютер связан с каждым, то есть имеем полностью связанную систему. Тогда процесс вычислений организуется согласно рис. 2.

Кіі

К21

Ки

Кі2

К22

К;2

Кіз

К23

Кз

Кш

KsN

Рис. 1. Схема параллельной вычислительной системы Шаг 1. Вычислить Ку и 3^ по формулам

Зг,! = 32,1 Ку ; 8г,1= ,7=1, 2, N

и переслать Ку нац всем Л/ — 1 компьютеру.

Шаг 2. Вычислить Кг,/ и аз,; по формулам

с>з,./ = аз,1 Ку +аз,2 Кг,/; <%,у = у^ +Лй,у,7=1, 2, ..., А/

и переслать Кг,/ и ^з,у всем (А/-1) компьютерам.

Шаг в. Вычислить К;,/.

Шаг (в+1). Вычислить ^и+1) = Ь1К1,у + ЬгКг.у + ... + Ь0К0,;, 7 = 1,2, ..., А/. Шаг (в+2). Сохранить = ^и+1).

Шаг (в+З). Переслать _у^и+1)всем (Л/-1) компьютерам.

Рис. 2. Вычисление коэффициентов и вектора приближенного решения при р = N

На последнем шаге в каждом компьютере вычисляется и сохраняется 7-я компонента вектора /й+1). Число пересылок для расчета в одном узле сетки \/\/п составит - 1)2« 0(вN2). При высокой размерности исходной задачи (1) затраты на пересылку являются определяющими для времени решения (1). Сокращение затрат в общем случае можно уменьшить только укрупнением схемы параллельных вычислений, то есть за счет увеличения числа коэффициентов, вычисляемых на одном компьютере. Рассмотрим в качестве примера схему, которая показана на рис.3. Здесь для определенности размерность N кратна числу р, то есть N = кр. В каждом компьютере размещено и вычисляется последовательно по к компонент векторов коэффициентов, Кі(п), К2(п),., ^(п\ Вычисления реализуются по правилу, показанному на рис. 4.

Кіі, Кі2, ..., Кік К>і, К22, ..., К2к

Ks1, Ks2, ■■■, Ksк

Кі,к+і, ..., К|,2к К2,к+і, ..., К2,2к

Кз,к^і, ■ ■■, ^,2к

Кі, (р-і)к+і, , ..., Кі,рк Кг, (р-і)к+і, , ..., К2,рк

К

■в, (р-і)к+і,

К

з, рк

Рис. 3. Укрупненная схема параллельной вычислительной системы Шаг 1. Вычислить К^ и б^ по формулам

б<и) - у (,г) + ha К™ z =1,2,..., р\/г - (=■-\)к +1,2,..., zk

’ Z Z ’ ’ Z

и переслать К^ и бвсем (р-1) компьютерам.

’ Z ’ z

Шаг 2. Вычислить Kи б^} по формулам

б!,”,’ =У, {"'+1г^а,n JKf] ~(z-l)k + \,2,...,zk

’ z z J’ z

j —1

и переслать K(nJ и б(nJ всем (р-1) компьютерам.

Шаг s. Вычислить K

(и)

, l 9 z

Шаг (s+1). Вычислить

yL^=yM + '£b¡K^,:=U2,...,piL=(:-m + U2,...,:k.

. - - _ I

z z

m-1

Шаг (s+2). Сохранить = ÿ)n+l).

Z Z

Шаг (s+3). Переслать j)/(”+1)BceM (р-1) компьютерам.

z

Л

Рис. 4. Вычисление коэффициентов и вектора приближенного решения при N = kp

i i Г,(,г + 1) т

На последнем шаге в каждом компьютере вычисляется и сохраняется своя часть вектора у .Таким образом, после параллельного вычисления коэффициентов параллельно определяется такое же количество компонент вектора приближенного решения у{п+1) ■ Число пересылок для одного узла сетки wn составит s(p-1)2«0(sp2). Заметим, что в отдельных случаях, зависящих от типа и вида исходной системы (1), может быть достигнута максимальная степень параллелизма и соответственно сокращение временных затрат выполнения вычислений по разработанному алгоритму.

Заключение. Рассмотренные параллельные схемы явных методов типа Рунге-Кутта ориентированы на реализацию в многопроцессорных вычислительных системах кластерной архитектуры с использованием технологии MPI. MPI имеет в составе коммуникационные операции, попарные и коллективные обмены, средства организации виртуальных топологий. Как показали исследования представленных параллельных схем, для их реализация наиболее подходящими могут быть такие топологии, как кольцо, линейка, решетка и гиперкуб.

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

Литература

1. Воеводин В.В., Воеводин Вл.В. Параллельные вычисления. - СПб.: БХВ - Петербург, 2002. - 806 с.

2. Новиков Е.А. Явные методы для жестких систем. - Новосибирск: Наука, 1997. - 197с.

3. Хайрер Э., Ваннер Г. Решение обыкновенных дифференциальных уравнений. Жесткие и дифференциально-алгебраические задачи. - М.: Мир, 1999. - 685 с.

4. Jackson K., Norsett S. The potential for parallelism in Runge -Kutta methods. Part I: RK formulas in standart form // SIAM J. Numer. Anal. - 1996. - Vol. 32. - P. 49-82.

5. Hendrickson B., Kolda Tamara G. Graph partitioning models for parallel computing // Parallel Computing. -2002. - Т. 26. - № 12. - P. 181-197.