CC BY
11
www.volsu.ru
см
«
га
н <
03
CJ
о К
к
со
и
о а С
03 03
га
s
tr
001: https://doi.Org/10.15688/mpcm.jvolsu.2021.2.4
УДК 519.6 Дата поступления статьи: 08.02.2021
ББК 26.22 Дата принятия статьи: 08.04.2021
ЛОКАЛЬНО-ДВУМЕРНЫЕ СХЕМЫ РАСЩЕПЛЕНИЯ ДЛЯ ПАРАЛЛЕЛЬНОГО РЕШЕНИЯ ТРЕХМЕРНОЙ ЗАДАЧИ ТРАНСПОРТА ВЗВЕШЕННОГО ВЕЩЕСТВА1
Александр Иванович Сухинов
§ Член-корреспондент РАН, доктор физико-математических наук, профессор,
заведующий кафедрой математики и информатики, Донской государственный технический университет к sukhinov@gmail.com
https://orcid.org/0000-0002-5875-1523
пл. Гагарина, 1, 344003 г. Ростов-на-Дону, Российская Федерация
Александр Евгеньевич Чистяков
Доктор физико-математических наук, профессор кафедры математики и информатики,
к Донской государственный технический университет
https://orcid.org/0000-0002-8323-6005
Ш cheese_05@mail.ru
о к
О пл. Гагарина, 1, 344003 г. Ростов-на-Дону, Российская Федерация
Валентина Владимировна Сидорякина
X Кандидат физико-математических наук, доцент, заведующая кафедрой математики,
< Таганрогский институт имени А.П. Чехова (филиал)
Ростовского государственного экономического университета (РИНХ)
g cvv9@mail.ru
https://orcid.org/0000-0001-7744-015X @ ул. Инициативная, 48, 347936 г. Таганрог, Российская Федерация
Софья Владимировна Проценко
Аспирант кафедры математики и информатики, Донской государственный технический университет rab55555@rambler.ru
пл. Гагарина, 1, 344003 г. Ростов-на-Дону, Российская Федерация
Ася Михайловна Атаян
Аспирант кафедры математики и информатики, Донской государственный технический университет atayan24@mail.ru
пл. Гагарина, 1, 344003 г. Ростов-на-Дону, Российская Федерация
Аннотация. Рассматривается 30-модель транспорта взвесей в прибрежных морских системах, учитывающая множество факторов, среди которых гидравлическая крупность или скорость осаждения частиц, распространение взвеси, оседание, интенсивность распределения источников взвешенного вещества и др. Разностные операторы диффузионного переноса в горизонтальных и вертикальном направлениях для данной задачи обладают существенно различными характерными пространственно-временными масштабами процессов, а также спектрами. При типичной дискретизации, применительно к мелководным системам Юга России (Азовское море, Цимлянское водохранилище), шаги по горизонтальным направлениям составляют 200-1000 м, коэффициенты турбулентного обмена (турбулентной диффузии) (103-104)м2/с; по вертикальному направлению — шаги 0,1 м-1 м, а коэффициенты микротурбулентного обмена по вертикали — (0,1-1)м2/с. Если ориентироваться на использование явных локально-двумерных — локально-одномерных схем расщепления, то допустимые значения шага по времени для двумерной задачи составят порядка 10-100 с, а для одномерной задачи по вертикальному направлению — 0,1-1 с. Сказанное мотивирует к построению аддитивной локально-двумерной — локально-одномерной схемы расщепления по геометрическим направлениям. В работе приведено описание параллельного алгоритма, использующего для аппроксимации двумерной задачи диффузии-конвекции по горизонтальным направлениям и одномерной задачи диффузии-конвекции по вертикальному направлению как явными, так и неявными схемами. Двумерная неявная задача диффузии-конвекции по горизонтальным направлениям численно решается адаптивным попеременно-треугольным методом. Численная реализация одномерной задачи диффузии-конвекции по вертикальному направлению осуществляется последовательным методом прогонки для серии независимых на данном слое одномерных трехточечных задач по вертикальному направлению. Для повышения эффективности параллельных расчетов также выполнена декомпозиция расчетной пространственной сетки и всех сеточных данных по одному или двум пространственным направлениям — в горизонтальных направлениях. Проведено сравнение полученных алгоритмов с учетом допустимых величин шагов по времени и реальных временных затрат на выполнение вычислений и обменов информацией на каждом временном слое.
Ключевые слова: транспорт взвешенного вещества, численные методы, схемы расщепления, локально-двумерные схемы расщепления, параллельные алгоритмы.
Введение
В статье рассматривается 30-модель диффузии-конвекции-осаждения взвеси в водной среде и ее параллельная численная реализация на основе явных и неявных схем, а также расщепления по геометрическим направлениям с учетом особенностей прибрежной системы.
Реальная трехмерная турбулентность в прибрежных системах обладает рядом специфических особенностей, описанных, например, в работах [4; 9; 12]. Строго говоря, она может считаться локально-однородной и изотропной лишь в масштабах, значительно меньших глубины акватории, так как стратификацию воды хорошо перемешанного мелкого водоема ограничивает максимальный вертикальный размер турбулентных вихрей. Считается, что по последней причине процесс турбулентного перемешивания можно разделить на горизонтальный турбулентный обмен и вертикальную турбулентную диффузию. Операторы конвективного и диффузионного переноса в горизонтальных и вертикальном направлениях для задач транспорта взвеси обладают существенно различными физическими и спектральными свойствами [7; 11; 18]. В связи со сказанным, для построения экономичных алгоритмов оперативного прогноза распространения взвеси (загрязнения водоема) целесообразно применять локально-двумерные схемы (ЛДС) расщепления [1; 16]. Ранее было показано для пространственно-трехмерных начально-краевых задач для уравнений параболического типа, что применение ЛДС позволяет получать малозатратные по количеству арифметических операций алгоритмы решения параболических сеточных уравнений, которые являются экономичными по числу операций и времени, необходимого для осуществления обменов между процессорами [17].
1. Непрерывная 3D-модель диффузии-конвекции-осаждения взвесей
Будем использовать прямоугольную декартову систему координат Охуг, где ось Ох проходит по поверхности невозмущенной водной поверхности и направлена в сторону моря. Пусть К = Н + п — общая глубина акватории, [м]; Н — глубина при невозмущенной поверхности водоема, [м]; п — возвышение свободной поверхности относительно геоида, [м].
В модели осуществим переход от г-координатной системы к 0-координатной системе, для которой используем декартову систему координат в горизонтальной плоскости, а в качестве вертикальной переменной — безразмерную переменную 0, 0 € [0; 1].
В 0-координатной системе водный столб разделен на одинаковое количество слоев в каждой точке независимо от глубины, поэтому при использовании «новой» координатной системы решаются некоторые проблемы, связанные с добавлением и вычитанием слоев.
При переходе к 0-координатной системе используем формулу:
(а - Ъ)(г - п) (1)
0 = а----, Х0 = х, у0 = у, (1)
К
где 0 = а = 0 на свободной поверхности водоема (верхняя граница), 0 = 6 = 0 на дне (см. рис. 1).
Рис. 1. 0-координатная система
Далее вместо выражений (1) будем использовать соотношения
* - П
0
h
Х0
У0 = у.
(2)
Пусть в водном объеме V = {0 < х < Ьх, 0 < у < Ьу, 0 < 0 < 1} находятся частицы взвеси, которые в точке (х,у, 0) и в момент времени Ь имеют концентрацию с = с(х,у, 0, ¿), [мг/л]; Ь — временная переменная, ¿е [0; Т] [сек].
Уравнение, описывающее поведение частиц, будет выглядеть следующим образом [6; 13]:
дс д(ис) d(vc) a-bd((w + wg) с)
dt + дх
д
h
д0
(д2с д2с\ /а -Ь\2д ( дс\ ^ _
чдх^+W*)Л—)д00{уß00)+R (3)
Здесь и, v, w — компоненты вектора U скорости движения жидкости, [м/с]; wg — гидравлическая крупность или скорость осаждения частиц, [м/с]; ц, v — коэффициенты горизонтальной и вертикальной диффузии частиц, [м2/с]; F = F(х,у, 0, t) — функция, описывающая интенсивность распределения источников взвешенного вещества.
В совокупности с граничными условиями для функции концентрации частиц решение уравнения (3) позволяет определять потоки взвешенного вещества как по направлению к берегу, так и вдоль берега.
В качестве водного объема V рассматриваем трехмерную цилиндрическую область, верхнее основание которой лежит на свободной поверхности, а нижнее основание является частью донной поверхности. Отметим, что постоянная высота области V обусловлена переходом от координаты z по вертикали к координате 0, 0 G [0; 1]. Пусть Г = S U Stop U Sbottom, — граница области V; S — боковая граничная поверхность; Stop, Sbottom — части свободной поверхности и поверхности дна соответственно. Боковая поверхность S области V образуется в результате движения образующей по замкнутой линии, являющейся контуром аппроксимирующим береговую линию или «жидкую границу» рассматриваемой области. Областью задания уравнения (3) считаем
QT = V х {0 <t <Т}, V = V U Г.
Добавим к уравнению (3) начальные и граничные условия (предполагая, что осаждение частиц на дно необратимо):
• начальные условия при £ = 0
с(х,у, 0, 0) = Со(х,у, 0); (4)
• граничные условия на боковой границе $ в любой момент времени Ь
^ = 0, если (иг,п) < 0, (5)
Оп V /
дс иг тт _
--с, если иг,п) >0
(Ur,n) > 0, (6)
дп ц
где п — внешняя нормаль к границе области S, Uг — вектор скорости движения жидкости на границе S, иг — проекция вектора скорости Uг на направление нормали п на границе области S;
• граничные условия на поверхности воды Stop
I = 0; (7)
• граничные условия на дне Sbottom
| = -W С (8)
00 V
считая, что донная поверхность не «сильно» искривлена.
В работах [3;19] исследованы условия корректности задачи транспорта взвесей, в случае многокомпонентного гранулометрического состава частиц, а потому и задачи (3)-(8), как ее частного случая при условиях гладкости функции решения
о{х,у, 0,t) е C2(QT) ПС (QT), grad с еС (QT) и необходимой гладкости границы области.
2. Метод построения аддитивной локально-двумерной схемы расщепления
Пусть функция с имеет ограниченные производные до четвертого порядка включительно по декартовым переменным х, у, 0 и до второго порядка включительно по временной переменной t, а функции F, F\, F2, для которых F = F\ + F2, имеют ограниченные производные до второго порядка включительно по х, , 0 и первого порядка по t в V.
Обозначим и1 = и, и2 = v, и3 = w + wg.
Представим пространственно-трехмерные операторы конвективного переноса в форме:
^ д(и\с) д(и2с) а — b д(и3с) с дх ду h д0
1 / д (щ с)
2 V дх
\ 1 /д(щ с) дс\ 1а -bf
+ щаХ) + + щду) + 2~h~ \
( д2с д2с \ (а - Ь\2д (
DC = Ц удхХ2 + ду2) Л~)дв К
1а -bf ди3 дс
2V~6Q дв
)
дв
(9) (10)
где Ьс — линейный кососимметрический оператор. Вводя следующие обозначения
Lc12C "2V дх
д
х
D\2C = Ц
1 / д(щс)
П ду д2 д2
д
+ и2 ду), Lc3°
i д^с д^с ^ _/ а—b \2д_ ( дс \
Vöx2 + ду2) , 3С v,
а — b f ди3
д
2h I дв +Щздв
ч дх2 ду2
для формул (9), (10) можем записать равенства
Ьсс = Ьс12С + Ьс3 с, Вс = Б12С + Дз с. Используя введенные операторы Ьсс, Дс, уравнение (3) представим в виде
дс
— = -Ьсс + Дс + Д, ¿е [0; Т].
(11) (12)
(13)
Будем для простоты рассматривать случай, когда функция источника действует на поверхности водоема, то есть выполнено условие Д(х,у, 0, ¿) = Д(х,у, ¿).
На временном отрезке [0; Т] построим равномерную сетку шт с шагом т, то есть множество точек вида шт = {= пт, п = 0, М, Мт = Т}.
На введенной временной сетке шт заменим трехмерную задачу (13), (4)-(8) цепочкой двумерных-одномерных задач, связанных по начальным и финальным данным — значениям сеточных функций С(\) и С(2):
д
(1)
д
- LC12C(1) + Di2C(i) + Fl, tn < t< tn+i, n = 0,N - 1,
С(1)(х,У, 0, 0) = Co(х,y, 0), С(2)(х,у, 0, tn) = С(1)(х,у, 0, tn+i), n = 0,N - 1
(14)
(15)
(16)
(2)
д
-LC3C(2) + D3C(2) + F2, tn < t < tn+i,
n
0, N 1,
С(1)(х,у, 0, tn) = С(2)(х,у, 0, tn),
n
1,N - 1.
(17)
(18)
Уравнения (14)-(16) и (17), (18) дополняются граничными условиями (5), (6) на боковой границе и граничными условиями (7), (8) на поверхностях воды и дна соответственно.
Решением этой задачи является функция с(х,у, 0, Ьп) = С(2)(х,у, 0, Ьп), п = 1,М. В области V построим равномерные сетки
whi2 = {хг = гhx, yj = jhy;i = 0,Nx,j = 0,Ny; Nxhx = Lx,Nyhy = Ly}
и
шНэ = [0к = к!10; к = 0, N0; Щк0 = 1},
где Нх, ку, к0 — шаги по пространству; N, N — количество узлов по пространственным координатам х, у, 0; и образуем пространственную сетку ш^, = ш^12 х
Дискретизация задач (14)-(28), входящих в расщепленную цепочку, осуществляется на пространственно-временной сетке = х шт с помощью интегро-интерполяционного метода. Выражения Ьс12С(1), Ьс3С(2), ^12^(1), ^эС(2) аппроксимируются соответственно разностными выражениями К12с, Кэс, Л12с, Лэс интегро-интерполяционным методом на прямоугольной равномерной сетке с использованием функции заполненности ячеек. Имеем
С(1)(х,у, 0, 0) = со(х,у, 0), (х,у, 0) е шн, при п = 0, (19)
С(\) = cf2), (х,у, 0, tn) Е whr, п = 1,N - 1, (20)
га+1 _ -n
ГШ_ЧЦ-Л .-.n+1 ^ _™+1
т
'""(2) = С(1)
с?*1 - cn
Л^!1 -К^с?*1 + ф?, (х,у, 0) Е œhT, п = 0, N - 1, (21)
с-2) = с?*1, (х,у, 0) G WhT, п = 0, N - 1, (22)
(2) - (2)
- Л ^n+1 ТУ П+1
т
Лэс?**1 -КэС?**1 + фП, tn < t< tn+1, п = 0, N. (23)
Уравнения (22)-(24) дополняются выражениями с соответствующей аппроксимацией граничных условий (5)-(6), а уравнения (22)-(23) — условий (7)-(8).
Для решения трехточечных разностных задач (25)-(26) применяются последовательные алгоритмы прогонки. Таким образом, на каждом временном шаге последовательно решаются две начально-краевые разностные задачи (19)-(21) с соответствующими граничными условиями на боковой поверхности и (22)-(23) с соответствующими граничными условиями на свободной поверхности и на дне. В качестве методов решения двумерной задачи на верхнем переменном слое будем использовать адаптивный попеременно-треугольный метод, который построен для случая ограниченного значения сеточного числа Пекле — Ре^, Ре^, < 2; для морских и прибрежных систем данное ограничение выполняется с достаточно большим запасом.
3. Метод решения сеточных уравнений
Наиболее ресурсоемкой частью программной реализации является численная реализация системы двумерных сеточных уравнений при аппроксимации ЛДС неявными схемами (схемами с весами) или при полной аппроксимации неявными схемами решения системы линейных алгебраических уравнений (СЛАУ). Для ускорения вычислений был разработан параллельный алгоритм, основанный на адаптивном попеременно-треугольном методе (ПТМ) решения СЛАУ и использующий технологию MPI для осуществления информационного обмена между вычислителями.
В конечномерном гильбертовом пространстве H рассматривается задача об отыскании решения операторного уравнения
Ах = f, А : H ^ H, (24)
где А — линейный, в общем случае несамосопряженный, положительно определенный оператор (А > 0). Для нахождения решения задачи (24) будем использовать неявный итерационный процесс
тт+1 _ —т
В--— + Ахт = f. (25)
Цт+1
В уравнении (25) m — номер итерации, Пт+i > 0 — итерационный параметр, а В — некоторый обратимый оператор, который называется предобуславливателем или стабилизатором. Обращение оператора В в (25) должно быть существенно проще, чем непосредственное обращение исходного оператора А в (24). Оператор В построен, исходя из аддитивного представления оператора А0 — симметричной части оператора А и его кососимметричной части
Ао = Ri + R2, Ri = R*, (26)
где А = А0 + А1, А0 = А*, А1 = —А1 на основе симметричной части А0 оператора А.
Оператор-предобуславливатель запишется в следующем виде
В = (D + wR1)D-1(D + wR2), D = D*, ш > 0, (27)
где D — некоторый оператор. Соотношения (26), (27) задают адаптивный ПТМ решения сеточных уравнений [14], если определены операторы R1, R2 и указаны способы определения параметров цт+1, ш и оператора D.
4. Параллельная реализация
Программная реализация модели гидродинамики была осуществлена на языке C++ с использованием технологии MPI, что позволило провести ряд вычислительных экспериментов на многопроцессорной вычислительной системе Научно-технологического университета «Сириус». Использовалась часть кластера, построенная на основе открытого стека приложений OpenHPC, содержащая 1440 процессорных ядер по 2,3 ГГц и 10,24 Тб оперативной памяти.
При параллельной реализации использованы методы декомпозиции сеточных областей для вычислительно трудоемких задач диффузии-конвекции, учитывающие архитектуру и параметры многопроцессорной вычислительной системы [5; 8; 10; 15; 20]. Декомпозиция расчетной двумерной области выполнена по двум пространственным переменным х, у (см. рис. 2), а также использовалась декомпозиция по одному пространственному направлению (одной горизонтальной координате).
Рис. 2. Декомпозиция двумерной сеточной области
На первом шаге вычислений первый процессор обрабатывает верхний слой. Затем осуществляется передача перекрывающихся элементов смежным процессорам. На
следующем шаге первый процессор обрабатывает второй слой, а его соседи — первый. Передача элементов после расчета двух слоев первым процессором показана на рисунке 3.
Рис. 3. Схема для расчета вектора решений для системы с нижней треугольной матрицей
В схеме для расчета вектора ут только первый процессор не требует дополнительной информации и может независимо от других процессоров вести обработку своей части области, остальные процессоры ждут результатов от предыдущего процессора, пока тот не передаст вычисленные значения сеточных функций для узлов сетки, располагающихся в предшествующих позициях данной строки. Процесс продолжается до тех пор, пока не будут рассчитаны все слои. Далее вычисляются скалярные произведения
(шт шт) = ^(шт)2, (28)
г
(Ашт, шт) = (29)
г
(Ашт, Ашт) = ^(Ашт)г2. (30)
г
Каждый процессор считает свою часть скалярных произведений (28)-(30), затем осуществляется передача данных от всех процессоров (ядер) ко всем. Значение скалярных произведений будет равно сумме всех частичных сумм. С помощью формулы (28) осуществляется переход на следующую итерацию. Представим результаты численных расчетов.
Были использованы расчетные сетки nN'х х kNу х N0 , где параметры п и к определяют способ декомпозиции по двум (одному) из координатных горизонтальных направлений, п^х = кNу = N, пк = р.
В таблице 1 приведены зависимости времени моделирования от распределения процессоров по пространственным координатам (х,у, 0). Для расчетных сеток с числом узлов по трем координатным направлениям 200x200x100 и 1000x1000x100 узлов при разбиении расчетной области на р частей по одному пространственному направлению.
В таблице 2 представлены полученные значения ускорения и эффективности в зависимости от количества процессоров для параллельного варианта попеременно треугольного метода. Здесь использованы обозначения: р — количество процессоров; п, к — количество подобластей по пространственным координатам х и у соответственно.
Таблица 1
Результаты работы алгоритма — время решения в секундах — с использованием технологии MPI для различного числа вычислителей на двух расчетных сетках(схема расщепления — неявные ЛДС — ЛОС)
Сетка/число процессоров 1 2 4 8 12 16 20 24
200x200x100 1,3167 0,7747 0,3401 0,2053 0,1392 0,1352 0,1118 0,0974
1000x1000x100 55,303 30,037 16,640 6,8903 4,6011 4,421 3,268 2,9
Таблица 2
Результаты работы параллельного алгоритма, построенного на основе декомпозиции расчетной области по двум пространственным направлениям (неявные схемы расщепления ЛДС — ЛОС)
Р п к 200x200x100 1000x1000x100
время ускорение эффективность время ускорение эффективность
16 16 1 0,1352 9,74 0,6757 4,421 12,51 0,7818
16 8 2 0,1232 10,69 0,6854 3,52 15,71 0,9821
16 4 4 0,1123 11,72 0,7881 3,499 15,8 0,9878
20 20 1 0,1118 11,78 0,5891 3,268 16,92 0,846
20 10 2 0,1104 11,93 0,5965 3,203 17,27 0,8633
20 5 4 0,09686 13,6 0,6798 3,153 17,54 0,877
24 24 1 0,09744 13,51 0,5631 2,9 19,07 0,7946
24 12 2 0,09134 14,42 0,6007 2,846 19,43 0,8096
24 8 3 0,08365 15,74 0,6559 2,831 19,53 0,8139
24 6 4 0,07407 17,78 0,7407 2,77 19,97 0,8319
Коэффициент ускорения на прямоугольной сетке 1000x1000x100 узлов при числе использованных процессоров, равном 24, составил 19,07 раза (время исполнения для одного процессора — 55,30 с., для 24 — 2,9 с). В таблице 3 представлены зависимости времени выполнения одного шага алгоритма ПТМ от числа процессоров для двух расчетных сеток, содержащих: 500x500x100 и 1000x1000x60 узлов при декомпозиции расчетной области по одному пространственному направлению.
Таблица 3
Время исполнения алгоритма в с. с использованием технологии MPI для различного числа вычислителей для двух расчетных сеток (ПТМ)
Сетка 1 2 4 8 12 16 20 24
500x500x100 15,652 8,1459 4,1752 2,2188 1,5997 1,3111 1,0814 0,9720
1000x1000x100 38,16 19,690 10,210 5,3694 3,8203 3,1618 2,6823 2,5091
Результаты расчета ускорения и эффективности в зависимости от количества процессоров для параллельного варианта попеременно треугольного метода приведены в таблице 4.
Коэффициент ускорения алгоритма при использовании прямоугольной сетки, содержащей 1000x1000x60 узлов и числе процессоров, равном 24, составило 15,2 (время исполнения для одного процессора — 38,16 с., для 24 — 2,51 с). Подводя итоги, можно сказать, что схема расщепления, использующая явную ЛДС и неявную ЛОС имеет преимущество при использовании относительно грубых сеток по горизонтальным направлениям для моделирования процесса переноса взвесей в прибрежных системах (500-1000 м, характерное число узлов по каждому из горизонтальных координатных направлений 200-400 м), по вертикальному направлению порядка 100 м. При использовании более подробных сеток (характерный шаг сетки 50-300 м по каждому из горизонтальных координатных направлений, характерное число узлов по каждому из горизонтальных координатных направлений — 1000 м и более) наилучшими являются неявные схемы расщепления ЛДС — ЛОС. Данные выводы применимы для относительно небольшого числа процессоров (немногие десятки) и в существенной мере зависят от архитектурных особенностей конкретной вычислительной системы.
Таблица 4
Время выполнения параллельного алгоритма адаптивного ПТМ, построенного на основе декомпозиции расчетной области по двум пространственным направлениям, для численной реализации двумерной разностной схемы, входящей в аддитивную схему ЛДС — ЛОС
р п к 500x500x100 1000x1000x60
время ускорение эффективность время ускорение эффективность
16 16 1 1,311 11,94 0,7461 3,162 12,07 0,7543
16 8 2 1,288 12,15 0,7594 3,004 12,7 0,794
16 4 4 1,284 12,19 0,7619 2,987 12,78 0,7985
20 20 1 1,081 14,47 0,7237 2,682 14,23 0,7113
20 10 2 1,033 15,16 0,758 2,547 14,98 0,749
20 5 4 1,023 15,3 0,7652 2,499 15,27 0,7635
24 24 1 0,9721 16,1 0,6709 2,509 15,21 0,6337
24 12 2 0,9629 16,26 0,6773 2,246 16,99 0,708
24 8 3 0,9462 16,54 0,6893 2,242 17,02 0,7091
24 6 4 0,9316 16,8 0,7001 2,241 17,03 0,7095
Заключение
В статье рассмотрено построение и исследование параллельных алгоритмов для реализации 30-модели транспорта взвесей в прибрежных морских системах на основе ЛДС расщепления. Проведено сравнение параллельных алгоритмов, использующих ЛДС с явной и неявной аппроксимацией двумерных задач, входящих в аддитивную схему. Для численной реализации применен адаптивный ПТМ, разработанный ранее в авторском коллективе для разностных задач уравнений с несамосопряженным операто-
ром для ограниченного числа Пекле (Pe^ < 2). Для моделирования динамики транспорта взвеси прибрежной системы разработан программный комплекс. На многопроцессорной вычислительной системе НТУ «Сириус» исследована параллельная реализация программного комплекса, выполненная на языке C++ с использованием технологии MPI. Результаты работы параллельного алгоритма, построенного на основе декомпозиции расчетной области по двум пространственным направлениям, позволили оценить его эффективность.
ПРИМЕЧАНИЕ
1 Исследование выполнено при финансовой поддержке РФФИ в рамках научного проекта № 20-01-00421.
СПИСОК ЛИТЕРАТУРЫ
1. Сухинов, А. И. Двумерные схемы расщепления и некоторые их приложения / А. И. Сухинов. - М. : Макс Пресс, 2005. - 408 с.
2. Сухинов, А. И. Математическое моделирование транспорта наносов в прибрежных водных системах на многопроцессорной вычислительной системе / А. И. Сухинов, А. Е. Чистяков, Е. А. Проценко // Выч. мет. программирование. — 2014. — Т. 15, вып. 4. — C. 610-620.
3. Сухинов, А. И. Построение и исследование корректности математической модели транспорта и осаждения взвесей с учетом изменения рельефа дна / А. И. Сухинов, В. В. Си-дорякина // Вестник Донского государственного технического университета. — 2018. — Т. 18, № 4. — C. 350-361. — DOI: https://doi.org/10.23947/1992-5980-2018-18-4-350-361.
4. Численная модель динамики поверхностных вод в русле Волги: оценка коэффициента шероховатости / А. В. Писарев, С. С. Храпов, Е. О. Агафонникова, А. В. Хоперсков // Вестн. Удмурт. ун-та. Матем. Мех. Компьют. науки. — 2013. — № 1. — C. 114-130.
5. Belotserkovskii, O. M. Decomposition method applied to the solution of problems of viscous incompressible fluid dynamics / O. M. Belotserkovskii, V. A. Gushchin, V. V. Shchennikov // Computational Mathematics and Mathematical Physics. — 1975. — Vol. 15. — P. 197-207.
6. Coastal hydrodynamics in a windy lagoon / E. Alekseenko, B. Roux, A. Sukhinov, R. Kotarba, D. Fougere // Nonlinear Processes in Geophysics. — 2013. — Vol. 20. — P. 189-198. — DOI: https://doi.org/10.1016/j-.compfluid.2013.02.003.
7. Flood-based analysis of high-magnitude sediment transport using a non-parametric method / T. Francke, J. A. Lopez-Tarazon, D. Vericat, A. Bronstert, R. J. Batalla // Earth Surface Processes and Landforms. — 2008. — Vol. 33, iss. 13. — P. 2064-2077. — DOI: https://doi.org/10.1002/esp.1654.
8. Kovenya, V. M. Application of Splitting Algorithms in the Method of Finite Volumes for Numerical Solution of the Navier-Stokes Equations / V. M. Kovenya, P. V. Babintsev // J. Appl. Ind. Math. — 2018. — Vol. 12. — P. 479-491. — DOI: https://doi.org/10.1134/S1990478918030080.
9. Leontyev, I. O. Coastal Dynamics: Waves, Moving Streams / I. O. Leontyev. — Moscow : GEOS Publ., 2001. — 272 p.
10. Milyukova, O. Yu. Parallel approximate factorization method for solving discrete elliptic equations / O. Yu. Milyukova // Parallel Computing. — 2001. — Vol. 27, iss. 10. — P. 1365-1379. — DOI: https://doi.org/10.1016/S0167-8191(01)00092-8.
11. Ouda, M. Development of a new multiphase sediment transport model for free surface flows / M. Ouda, E. A. Toorman // International Journal of Multiphase Flow. — 2019. — Vol. 117, iss. 3. — P. 81-102.
12. Predictive modeling in sediment transportation across multiple spatial scales in the Jialing River Basin of China / X. Liu, S. Qi, Y. Huang, Y. Chen, P. Du // International Journal of Sediment Research. - 2015. - Vol. 30, iss. 3. - P. 250-255. - DOI: https://doi.org/10.1016/j-.ijsrc.2015.03.013.
13. Samarskii, A. A. Finite-difference approximations to the transport equation. II / A. A. Samarskii, P. N. Vabishchevich // Differential Equations. - 2000. - Vol. 36. -P. 1069-1077. - DOI: https://doi.org/10.1007/BF02754509.
14. Sukhinov, A. I. Adaptive modified alternating triangular iterative method for solving grid equations with a non-self-adjoint operator / A. I. Sukhinov, A. E. Chistyakov // Mathematical Models and Computer Simulation. - 2012. - Vol. 4, iss. 4. - P. 398-409.
15. Sukhinov, A. I. Increasing efficiency of alternating triangular method based on improved spectral estimates / A. I. Sukhinov, A. V. Shishenya // Mathematical Models and Computer Simulations. - 2013. - Vol. 5. - P. 257-265. - DOI: https://doi.org/10.1134/S2070048213030101.
16. Sukhinov, A. I. Locally two-dimensional schemes for the approximation of the three-dimensional heat equation in toroidal coordinates / A. I. Sukhinov, V. S. Vasil'ev // Russian Mathematics (Iz. VUZ). - 1996. - Vol. 40, iss. 3. - P. 58-67.
17. Sukhinov, A. I. Numerical realization of the three-dimensional model of hydrodynamics for shallow water basins on a high-performance system / A. I. Sukhinov, A. E. Chistyakov, E. V. Alekseenko // Mathematical Models and Computer Simulations. - 2011. - № 3. -Article ID: 562. - DOI: https://doi.org/10.1134/S2070048211050115.
18. Sukhinov, A. I. Reconstruction of 2001 Ecological Disaster in the Azov Sea on the Basis of Precise Hydrophysics Models. / A. I. Sukhinov, A. A. Sukhinov // Parallel Computational Fluid Dynamics. - New York : Elsevier Science Inc., 2005. -P. 231-238. - DOI: https://doi.org/10.1016/B978-044452024-1/50030-0.
19. Sukhinov, A. I. Uniqueness of solving the problem of transport and sedimentation of multicomponent suspensions in coastal systems / A. I. Sukhinov, A. A. Sukhinov, V. V. Sidoryakina // Journal of Physics: Conference Series. - 2020. - Vol. 1479. - Article ID: 012081. - DOI: https://doi.org/10.1088/1742-6596/1479/1Z012081.
20. Two-Phase Porous Media Flow Simulation on a Hybrid Cluster. / M. Trapeznikova, B. Chetverushkin, N. Churbanova, D. Morozov // Large-Scale Scientific Computing. LSSC 2011. Lecture Notes in Computer Science. - Berlin, Heidelberg : Springer, 2012. -Vol. 7116. - P. 646-653. - DOI: https://doi.org/10.1007/978-3-642-29843-1_74.
REFERENCES
1. Sukhinov A.I. Dvumernye skhemy rasshchepleniya i nekotorye ikh prilozheniya [Two-Dimensional Splitting Schemes and Some of Their Applications]. Moscow, Maks Press Publ., 2005. 408 p.
2. Sukhinov A.I., Chistyakov A.E., Protsenko E.A. Matematicheskoe modelirovanie transporta nanosov v pribrezhnykh vodnykh sistemakh na mnogoprotsessornoy vychislitelnoy sisteme [Sediment Transport Mathematical Modeling in a Coastal Zone Using Multiprocessor Computing Systems]. Vych. met. programmirovanie, 2014, vol. 15, iss. 4, pp. 610-620.
3. Sukhinov A.I., Sidoryakina V.V. Postroenie i issledovanie korrektnosti matematicheskoy modeli transporta i osazhdeniya vzvesey s uchetom izmeneniya relyefa dna [Development and Correctness Analysis of the Mathematical Model of Transport and Suspension Sedimentation Depending on Bottom Relief Variation]. Vestnik Donskogo gosudarstvennogo tekhnicheskogo universiteta, 2018, vol. 18, no. 4, pp. 350-361. DOI: https://doi.org/10.23947/1992-5980-2018-18-4-350-361.
4. Pisarev A.V., Khrapov S.S., Agafonnikova E.O., Khoperskov A.V. Chislennaya model dinamiki poverkhnostnykh vod v rusle Volgi: otsenka koeffitsienta sherokhovatosti [Numerical Model of Shallow Water Dynamics in the Channel of the Volga: Estimation of Roughness]. Vestn. Udmurt. un-ta. Matem. Mekh. Kompyut. nauki, 2013, no. 1, pp. 114-130.
5. Belotserkovskii O.M., Gushchin V.A., Shchennikov V.V. Decomposition Method Applied to the Solution of Problems of Viscous Incompressible Fluid Dynamics. Computational Mathematics and Mathematical Physics, 1975, vol. 15, pp. 197-207.
6. Alekseenko E., Roux B., Sukhinov A., Kotarba R., Fougere D. Coastal Hydrodynamics in a Windy Lagoon. Nonlinear Processes in Geophysics, 2013, vol. 20, pp. 189-198. DOI: https://doi.org/10.1016/j.compfluid.2013.02.003.
7. Francke T., Lopez-Tarazon J.A., Vericat D., Bronstert A., Batalla R.J. Flood-Based Analysis of High-Magnitude Sediment Transport Using a Non-Parametric Method. Earth Surface Processes and Landforms, 2008, vol. 33, iss. 13, pp. 2064-2077. DOI: https://doi.org/10.1002/esp.1654.
8. Kovenya V.M., Babintsev P.V. Application of Splitting Algorithms in the Method of Finite Volumes for Numerical Solution of the Navier-Stokes Equations. J. Appl. Ind. Math., 2018, vol. 12, pp. 479-491. DOI: https://doi.org/10.1134/S1990478918030080.
9. Leontyev I.O. Coastal Dynamics: Waves, Moving Streams. Moscow, GEOS Publ., 2001. 272 p.
10. Milyukova O.Yu. Parallel Approximate Factorization Method for Solving Discrete Elliptic Equations. Parallel Computing, 2001, vol. 27, iss. 10, pp. 1365-1379. DOI: https://doi.org/10.1016/S0167-8191(01)00092-8.
11. Ouda M., Toorman E.A. Development of a New Multiphase Sediment Transport Model for Free Surface Flows. International Journal of Multiphase Flow, 2019, vol. 117, iss. 3, pp. 81-102.
12. Liu X., Qi S., Huang Y., Chen Y., Du P. Predictive Modeling in Sediment Transportation Across Multiple Spatial Scales in the Jialing River Basin of China. International Journal of Sediment Research, 2015, vol. 30, iss. 3, pp. 250-255. DOI: https://doi.org/10.1016/jj.ijsrc.2015.03.013.
13. Samarskii A.A., Vabishchevich P.N. Finite-Difference Approximations to the Transport Equation. II. Differential Equations, 2000, vol. 36, pp. 1069-1077. DOI: https://doi.org/10.1007/BF02754509.
14. Sukhinov A.I., Chistyakov A.E. Adaptive Modified Alternating Triangular Iterative Method for Solving Grid Equations with a Non-Self-Adjoint Operator. Mathematical Models and Computer Simulation, 2012, vol. 4, iss. 4, pp. 398-409.
15. Sukhinov A.I., Shishenya A.V. Increasing Efficiency of Alternating Triangular Method Based on Improved Spectral Estimates. Mathematical Models and Computer Simulations, 2013, vol. 5, pp. 257-265. DOI: https://doi.org/10.1134/S2070048213030101.
16. Sukhinov A.I., Vasil'ev V.S. Locally Two-Dimensional Schemes for the Approximation of the Three-Dimensional Heat Equation in Toroidal Coordinates. Russian Mathematics (Iz. VUZ), 1996, vol. 40, iss. 3, pp. 58-67.
17. Sukhinov A.I., Chistyakov A.E., Alekseenko E.V. Numerical Realization of the Three-Dimensional Model of Hydrodynamics for Shallow Water Basins on a High-Performance System. Mathematical Models and Computer Simulations, 2011, no. 3, article ID: 562. DOI: https://doi.org/10.1134/S2070048211050115.
18. Sukhinov A.I., Sukhinov A.A. Reconstruction of 2001 Ecological Disaster in the Azov Sea on the Basis of Precise Hydrophysics Models. Parallel Computational Fluid Dynamics. New York, Elsevier Science Inc., 2005, pp. 231-238. DOI: https://doi.org/10.1016/B978-044452024-1/50030-0.
19. Sukhinov A.I., Sukhinov A.A., Sidoryakina V.V. Uniqueness of Solving the Problem of Transport and Sedimentation of Multicomponent Suspensions in Coastal Systems. Journal of Physics: Conference Series, 2020, vol. 1479, article ID: 012081. DOI: https://doi.org/10.1088/1742-6596/1479/1Z012081.
20. Trapeznikova M., Chetverushkin B., Churbanova N., Morozov D. Two-Phase Porous Media Flow Simulation on a Hybrid Cluster. Large-Scale Scientific Computing. LSSC 2011. Lecture Notes in Computer Science. Berlin, Heidelberg, Springer, 2012, vol. 7116, pp. 646-653. DOI: https://doi.org/10.1007/978-3-642-29843-1_74.
LOCAL TWO-DIMENSIONAL SPLITTING SCHEMES FOR 3D SUSPENDED MATTER TRANSPORT PROBLEM PARALLEL SOLUTION
Alexander I. Sukhinov
Corresponding Member of the RAS, Doctor of Technical Sciences, Professor,
Head of Department of Mathematics and Physics,
Don State Technical University
sukhinov@gmail.com
https://orcid.org/0000-0002-5875-1523
Gagarina Sq., 1, 344002 Rostov-on-Don, Russian Federation
Alexander E. Chistyakov
Doctor of Technical Sciences, Professor, Professor of Department of Mathematics and Physics, Don State Technical University cheese_05@mail.ru
https://orcid.org/0000-0002-8323-6005
Gagarina Sq., 1, 344002 Rostov-on-Don, Russian Federation
Valentina V. Sidoryakina
Doctor of Technical Sciences, Professor, Head of Department of Mathematics, Chekhov Taganrog Institute (branch) Rostov State University of Economics cvv9@mail.ru
https://orcid.org/0000-0001-7744-015X
Initsiativnaya St, 48, 347924l Taganrog, Russian Federation
Sofya V. Protsenko
Postgraduate Student, Department of Mathematics and Physics,
Don State Technical University
rab55555@rambler.ru
Gagarina Sq., 1, 344002 Rostov-on-Don, Russian Federation
Asya M. Atayan
Postgraduate Student, Department of Mathematics and Physics,
Don State Technical University
atayan24@mail.ru
Gagarina Sq., 1, 344002 Rostov-on-Don, Russian Federation
Abstract. A 3D-model of suspended matter transport in coastal marine systems, which takes into account many factors, including the hydraulic size or the rate of particle deposition, the suspended matter propagation, sedimentation, the suspended matter sources distribution intensity, etc is considered. The difference operators of diffusion transport in the horizontal and vertical directions have significantly different characteristic spatiotemporal scales of processes, as well as spectra. With typical discretization, in relation to shallow-water systems in the Southern Russia (Azov Sea, Tsimlyansk reservoir), the steps in horizontal
directions are 200-1000 meters, the turbulent exchange coefficients (turbulent diffusion) (103-104) m2/sec; steps in vertical direction 0.1-1 m, and the vertical microturbulent exchange coefficients (0.1-1) m2/sec. If focus on the using of explicit local 2D-1D splitting schemes, then the permissible values of the time step for 2D problem will be about 10-100 seconds, and for 1D problem in the vertical direction 0.1-1 sec. This motivates to construct an additive local 2D-1D splitting scheme in geometric directions. The paper describes parallel algorithm that uses both explicit and implicit schemes to approximate the 2D diffusion-convection problem in horizontal directions and 1D diffusion-convection problem in the vertical direction. The 2D implicit diffusion-convection problem in horizontal directions is numerically solved by the adaptive alternating-triangular method. The numerical implementation of 1D diffusion-convection problem in the vertical direction is carried out by sequential run-through method for series of independent 1D three-point problems in the vertical direction on given layer. To increase the efficiency of parallel calculations, the computational spatial grid and all grid data are also decomposed into one or two spatial directions in the horizontal directions. The obtained algorithms are compared taking into account the permissible time steps values and the actual time spent on performing calculations and exchanging information on each time layer.
Key words: suspended matter transport, numerical methods, splitting schemes, locally two-dimensional splitting schemes, parallel algorithms.
