Численное решение прямых и обратных задач диффузии-конвекции на многопроцессорных системах для прогноза и ретроспективного анализа водных экосистем Текст научной статьи по специальности «Математика»

19. Chetverushkin B.N., Davydov A.A., Shil'nikov E.V. Simulating flows of incompressible and weakly compressible fluids on multicore hybrid computer systems, Computational Mathematics and Mathematical Physics, 2010, Vol. 50, No. 12, pp. 2157-2165.

Статью рекомендовал к опубликованию д.т.н., профессор Я.Е. Ромм.

Сухинов Александр Иванович - Южный федеральный университет; e-mail: sukhinov@gmail.com; 347928, г. Таганрог, пер. Некрасовский, 44; тел.: 88634315638; кафедра математического обеспечения суперкомпьютеров; профессор.

Шишеня Александр Владимирович - ассистент.

Sukhinov Alexander Ivanovich - Southern Federal University; e-mail: sukhinov@gmail.com; 44, Nekrasovsky, Taganrog, 347928, Russia; phone: +78634315638; the department of mathematical software supercomputers; professor.

Shishenya Alexsandr Vladimirovich - assistant.

УДК 532.5.031

Д.В. Лапин, А.Е. Чистяков, А.А. Сухинов

ЧИСЛЕННОЕ РЕШЕНИЕ ПРЯМЫХ И ОБРАТНЫХ ЗАДАЧ ДИФФУЗИИ-КОНВЕКЦИИ НА МНОГОПРОЦЕССОРНЫХ СИСТЕМАХ

ДЛЯ ПРОГНОЗА И РЕТРОСПЕКТИВНОГО АНАЛИЗА ВОДНЫХ

ЭКОСИСТЕМ*

Рассмотрены двумерные обратные задачи диффузии-конвекции, необходимость оперативного решения которых возникает при ретроспективном анализе техногенных и природных экологических катастроф. Алгоритм решения обратных задач основан на методе квазиобращения и последующего итерационного уточнения начального условия. Численные методы основаны на использовании схем расщепления по пространственным направлениям с последующим решением полученных систем линейных уравнений методом прогонки. Приведено описание предложенных параллельных алгоритмов и теоретические оценки ускорения и эффективности для систем с общей и распределенной памятью, а также гибридных систем. Поля течений для модели диффузии-конвекции рассчитываются при помощи трёхмерной гидродинамической модели, учитывающей стоки рек, геометрию дна и береговой линии, силу Кориолиса и прочие факторы. Общей особенностью предложенных моделей является учёт «заполненности» ячеек при расчёте на прямоугольной сетке. Заполненность отражает процент заполнения ячейки жидкостью, и позволяет более точно аппроксимировать криволинейную границу расчётной области по сравнению с бинарной функцией принадлежности ячейки расчётной области. На основе разработанного комплекса программ проведён численный эксперимент по восстановлению исходного поля концентрации загрязняющих веществ.

Перенос взвеси; обратные задачи; гидродинамика; численный эксперимент.

*Работа выполнена при частичной поддержке проектов Программы № 43 фундаментальных исследований Президиума РАН по стратегическим направлениям развития науки «Фундаментальные проблемы математического моделирования».

D.V. Lapin, A.E. Chistyakov, A.A. Sukhinov

NUMERICAL SOLUTION OF CONVECTION-DIFFUSION DIRECT AND INVERSE PROBLEMS ON THE MULTIPROCESSOR SYSTEMS FOR THE FORECAST AND RETROSPECTIVE ANALYSIS OF AQUATIC ECOSYSTEMS

The article considers two-dimensional convection-diffusion inverse problems, for which operative solution is needed for retrospective analysis of technogenic and natural environmental disasters. The algorithm for solving inverse problems is based on the method of quasi-inversion and subsequent iterative refinement of the initial conditions. Numerical methods are based on splitting schemes on spatial variables with subsequent solving of linear systems of equations by Thomas algorithm. There are described proposed parallel algorithms and given theoretical estimates of speed-up and efficiency for systems with shared and distributed memory, and hybrid systems. The flow fields for the diffusion-convection model are computed by three-dimensional hy-drodynamic model taking into account river flows, the shape of bottom and shoreline, Coriolis force and other factors. A common feature of the proposed models is taking into consideration of "fullness" of cells at the calculation on a rectangular grid. Fullness reflects the percentage of a cell filled with liquid, and can more accurately approximate the curved boundary of the computational domain compared with the binary membership function of the cell to computational domain. On the basis of developed software it was made a numerical experiment to recover the original field of pollutant concentration.

Sediments transport; inverse problem; hydrodynamics; numerical experiment.

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

Обратная эволюционная задача транспорта вещества. Будем рассматривать некорректную эволюционную задачу с обратным временем, получаемую из соответствующей прямой задачи заменой t на — t (переходом к обратному времени)

с;_ udx — vdy = —/ l—^y )y', (!)

где U = {u, v} - компоненты вектора скорости; /и - горизонтальная составляющая коэффициента турбулентного обмена.

Уравнение (1) рассматривается при следующих граничных условиях:

^ y, t\х,у y_r=0 (2)

и начальных условиях:

c(X У, t) t=0 = co (X У). (3)

Для приближенного решения некорректной задачи (1)-(3) будем использовать следующее уравнение:

c; — uc'x — vc'y =—/ )Х —(/cy )y' —a(c:xx )X —a(c'yyy ^. (4)

Данный подход был предложен в работе [1]. Уравнение (4) может быть записано в следующем виде:

1 '

с[ - ис[ - УС'У = (^ (-с - ас1 \ | + (-с - ас"№ ^ ^ . (5)

Для удобства последующей дискретизации уравнение (5) запишем в виде следующей системы уравнений:

c t - uc x - vcy =

f

У =И "x).' i + "у)y ] ' (6)

y

Sx =-C -ac"xx , (7)

^ =-с -ас'Уу. (8)

Выбор оптимального значения параметра регуляризации а определялся серией численных экспериментов, с использованием последовательности ак = а0дк, q > 0 [2]. Для каждого значения а оценивалось в норме сеточного пространства (щ), отклонение финального решения цепочки задач - обратная-прямая задача и начального распределения для обратной задачи на сетке Щ=ЮхХ (Оу (см. ниже).

Дискретизация обратной задачи транспорта вещества. Расчетная область вписана в прямоугольник. Для численной реализации дискретной математической модели поставленной задачи вводится равномерная сетка а = щхахха:

Щ =к = "К ,0 < п < N., I{ = к(Ы( }, щ = {хг = Шх, 0 < г < N, 1Х = }, Щ = ^ = К ,0 < ] < N у, ¡у = hyNу },

где п, г, ] - индексы по временной координате и пространственным координатным направлениям Ох, Оу соответственно; К, К, К - шаги по временной координате и пространственным координатным направлениям Ох, Оу соответственно; N , Nх, N - количество узлов по временной координате и пространственным координатным направлениям Ох, Оу соответственно; ¡(, ¡ж, I -длина

расчетной области по временной координате и пространственным координатным направлениям Ох, Оу .

Дискретные аналоги операторов конвективного и диффузионного переноса в случае граничных условий в форме Неймана могут быть записаны в следующем виде:

' / \ сг+1,] - сг,] , / \ сг,] - сг-1,]

исх « ( * К] -+ ( * и--,

г с — с с — с

(К )х * \] М+1/2,] 2 " - (\] г\2 '-1,] ,

где ^ - коэффициенты заполненности контрольных областей [2].

Дискретный аналог уравнения (6) с учетом приведенных выше аппроксимаций запишется в следующем виде:

_ сп сп+а1 _ Сп+а1 Сп+а1 _ Сп+а1

2П '—У2и 2П

СП+а2 _ СП+а2 СП+а2 _ СП+а2

_ (Ъ ), ,*и+Ч2-Ы1^-М--(Ъ )

V,.

2ку —- 12 2ку (9)

п+а / \ п+а / \ п+а / \ п+а

^ ) 1 _(^) (^) _(^ )

= (Ъ )',- Н+1/2,- ------(Ъ ),- Н_1/2,--^-~ +

(\ п+ 02 / \ п+ а ( \п+&г ( \ п+ а

+ (4 , и 5у_('у)п+ _(Ъ ) Н ^)п+ _('у)п+_1

+ (Ъ - и-+112 к2 (44 - и,-_112 к2 ,

где сп+а = ос"+ + (1 _ а)сп , а е [0,1] - вес схемы.

Дискретные аналоги уравнений (7)-(8) запишутся в следующем виде:

(4о- - = _(4о- С,- _ - С'2 С''3 _ (Ъ- ^ _

V Пх Пх У

(4о )и ('у = _(4о\-Си- _ « (4з )и _ (Ъ _

V пу пу

(10)

(11)

Для численного решения уравнения (5) с соответствующими начальными и граничными условиями используется схема расщепления по геометрическим направлениям - локально-одномерная. При этом для четных шагов веса схемы равны О = 1, а2 = 0, для нечетных - ах = 0, а2 = 1.

Прямое использование разностной схемы для численного решения уравнения (5) связано с численным решением сеточного эллиптического уравнения четвертого порядка, что достаточно сложно с точки зрения вычислительной трудоемкости. Поэтому целесообразно выполнить переход к аддитивной схеме расщепления по пространственным переменным [3] вида

2

^ + + _ Л>„ = 0.

т р=1

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

Математическая модель гидродинамики. Входными данными сконструированной выше обратной эволюционной задачи распространения и переноса вещества является поле вектора скорости, что требует, в свою очередь, построения математической модели движения водной среды. Исходными уравнениями гидродинамики мелководных водоемов являются [5-7]:

- уравнения движения (Навье-Стокса):

1 ' I \ ' '

и' + иы'х + Уи'у+ wu'z=--р'х+(ии'х)х +(ци'у) +(уи[}2 + 2^(у V 008^),

иу'х + уу'у + wvZ =_ 1 р'у +(ну'х)Х +(ну'у)у + 2Ш (12)

1 ' / \' '

+ + vw'х + ww'х =--р' +(им'х\ +(н''у) +(vw'z^ + 20и соъв;

- уравнение неразрывности в случае переменной плотности запишется следующим образом:

P't + P )x +My = 0 • (13)

где U = {u, v, w} - компоненты вектора скорости; p - превышение давления над гидростатическим давлением невозмущенной жидкости; р - плотность; Q - угловая скорость вращения земли; 0 - угол между вектором угловой скорости и вертикалью; ц, v - горизонтальная и вертикальная составляющие коэффициента турбулентного обмена.

Система уравнений (12)-(13) рассматривается при следующих граничных условиях:

- на входе (устье рек Дон и Кубань):

u(x,y,z,t) = u(t), v(x,y,z,t) = v(t), p'n(x,y,z,t) = 0, u'n(x,y,z,t) = 0,

- боковая граница (берег и дно):

PvM(u')n (X У, z, t) = -Tx (t), PvM(v')„ (X У, z, t) = -Ty (t),

p'n(x,y,z,t) = 0, u'„(x,y,z,t) = 0,

- верхняя граница:

PM(u')„(xУ,zt) = ~Tx(t) PM(v')„(xУ,zt) = ~Ty(ti (14)

w(x У, t) = -®- p't /Pg, p„ (x y, t) = 0,

- на выходе (Керченский пролив):

p„ (x,y, z, t) = 0, u'n (x,y, z, t) = 0,

где с - интенсивность испарения жидкости, тх , Т - составляющие тангенциального напряжения (закон Ван-Дорна), р - плотность взвеси.

Составляющие тангенциального напряжения для свободной поверхности:

Tx = PaCp ( Щ K|W|, ТУ = PaCp ( Щ Щу1Щ , где w - вектор скорости ветра относительно воды, Pa - плотность атмосферы, C (•) - безразмерный коэффициент.

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

Tx = PCp (u I)u|u I, ТУ = PCp IVN •

Рассмотренная ниже аппроксимация позволяет на основании измеренных (или вычисленных) пульсаций скоростей получить коэффициент вертикального турбулентного обмена, неоднородный по глубине[9], в соответствии с формулой:

v=с,2 Л21

rôûY fôVY

У ôz

+

ôz j

(15)

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

J

V

Построение и исследование дискретной модели. Расчетная область вписана в параллелепипед. Для численной реализации дискретной математической модели поставленной задачи гидродинамики вводится равномерная сетка:

wh = = пт, х = 1кх, у] = А, % = Щ;

0 < п < N, 0 < i < N, 0 < 7 < N, 0 < k < Nz; Nт = Т, Мхкх = 1Х, Муку = 1у, М2к2 = 1г }

где т - шаг по времени; кх, к , к - шаги по пространству; N - количество временных слоев; Т - верхняя граница по временной координате; Мх, , Мг -количество узлов по пространственным координатам; 1х, \ , /г - размеры параллелепипеда в направлении осей Ох, Оу и Oz соответственно.

Для решения задачи гидродинамики использовался метод поправки к давлению [7, 8]. Вариант данного метода в случае переменной плотности примет вид: й-и ' ' '

т У —У т

лЬ — МУ

г ' г

+ иу'х+уу'у+4>у'2 У^ + +(УК)2 —20.118шв,

I ' I

(16)

I '

Р" +р +р =--+ "-— +

г хх г уу V гг 2

Т Т Т Т

й-й 1 „, У — У 1 „, М>~М> 1 „,

- =--Рх, - =--Ру, - =--Р2,

т р т р т р

где {и, V, м} - компоненты вектора скорости; {и, V, М}, - компоненты

полей вектора скорости на «новом» и промежуточном временных слоях соответственно, и = (и + и )/ 2; р и р - распределение плотности водной среды на новом

и предыдущем временных слоях соответственно.

При построении дискретных математических моделей гидродинамики также учитывалась «заполненность» [4] контрольных ячеек, что позволяет повысить реальную точность решения в случае сложной геометрии исследуемой области за счет улучшения аппроксимации границы.

Через 01 jk обозначена «заполненность» ячейки (/, 7,к). Степень «заполненности» ячейки определяется давлением столба жидкости внутри данной ячейки. Если среднее давление в узлах, которые относятся к вершинам рассматриваемой ячейки, больше давления столба жидкости внутри ячейки, то ячейка считается заполненной полностью (^ к = 1). В общем случае «заполненность» ячеек можно вычислить по следующей формуле:

р + р + р + р

1 1,^ 1 1-1,7,к ^ 11,7-1,к ^ 11-1,7-1,к п1\

',],к рк2 где Р = р + рг - давление.

Погрешность аппроксимации математической модели равна о(г +||Ъ||2), где

= V+ Ъ + И'2 . Доказано сохранение потока на дискретном уровне разработанной гидродинамической модели, а также отсутствие неконсервативных дис-сипативных слагаемых, полученных в результате дискретизации системы уравнений. Достаточное условие устойчивости и монотонности разработанной модели определяется на основе принципа максимума [10] при ограничениях на шаг по

пространственным координатам: Ъ < |2[/м|, Ъ < |2[/, Ъ < \2у/^, или Яв < 2N, где Яв = V ■ 11[ - числа Рейнольдса, I - характерный размер области, N = тах^, N, N2}

Дискретные аналоги системы уравнений (16) решаются адаптивным модифицированным попеременно-треугольным методом вариационного типа [11-13].

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

Так как аддитивно-усредненная схема выбрана и для прямой и для обратной задачи, то принцип параллельного выполнения одинаков для обеих задач, разница заключается лишь в том, что использование аддитивно-усредненной схемы для прямой задачи приводит к необходимости решать трехдиагональные СЛАУ, а для обратной задачи - пятидиагональные СЛАУ. Такие СЛАУ решаются методом, соответственно, трех- и пятидиагональной прогонки [3].

Время, необходимое для решения таких СЛАУ на системе с общей памятью [14], имеющей к вычислителей можно оценить, соответственно, как: (7NM + 7MN)■ ^ и (23NM + 23MN)■ га

к к ' где ? - время выполнения одной арифметической операции; N и M - размеры

вычислительной области.

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

Параллельная реализация для систем с распределенной памятью. Параллельная реализация для систем с распределенной памятью отличается необходимостью выполнять транспонирование матрицы на решающем поле, рис. 1.

Рис. 1. Распределение неизвестных по узлам, транспонирование

При этом между шагами расщепления по пространственным переменным необходимо передать NM ■ (1 — 1/р) чисел с плавающей запятой двойной точности,

т.е. 8NM ■ (1 — 1/ р) байт. Здесь р - количество вычислителей в системе с распределенной памятью. Этот объем данных передается за р операций передач. Таким образом, время, затрачиваемое на передачу этих данных, составляет р ■ ^ + 8NM -(1 — 1/р)■ Хъ. Здесь ^ - латентность (время необходимое для инициализации передачи), а ^ - время, необходимое на передачу одного байта информации.

Для перехода к следующему шагу по времени для прямой задачи затрачивается время ((7NM + 7MN)■ ^)/ р на арифметические операции и

2(р ■ ^ + 8NM ■ (1 — 1/ р)- 1Ъ) на операции обмена, для обратной задачи время обмена такое же, а время выполнения арифметических операций составляет ((23 - NM + 23 - MN) - ^ )/р. После выполнения всех (К) шагов по времени для каждой из задач выполняется расчет невязки, требующий времени 5NM - ^ для выполнения арифметических операций и ^ + 8р - на выполнения операций обмена.

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

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

5 =_60 - ь - га - К+5 - ь - га_

р 5 - Ь - 1а + К-(4 - р - Г, + 8Ь-(1 — 1/р )-^ + 60 - Ь - та/(к - р )) + + 8р - ^ '

где Ь = M - N - характерный размер задачи.

Тогда эффективность параллельного выполнения на гибридной вычислительной системе

£ =

р к - р

Масштабируемость предложенного параллельного алгоритма может выполняться до тех пор, пока тт (Ы, N) > к - р.

Подставим в полученные оценки данные имеющейся в ЮФУ гибридной вычислительной системы. Вычислительная система ЮФУ обладает следующими характеристиками: ^ = 1,82 • 10—6 с, ^ = 8,16 -10—13 с, ¿о = 1,09 -10—10 с.

Подставляя эти данные, построим график для теоретической оценки ускорения (рис. 2, сплошная кривая), который достаточно хорошо согласуется с экспериментально полученными данными ускорения для гибридной супервычислительной системы ЮФУ (пунктирная кривая).

Рис. 2. Зависимости ускорения от числа процессоров: гладкая кривая - теоретическая зависимость; ломаная - экспериментальная

Результаты численных экспериментов. Численные эксперименты выполнялись на супервычислительной системе ЮФУ, представляющей собой гибридную вычислительную систему, состоящую из 128-ми узлов, объединенных сетью Infiniband 4x DDR ConnectX. Каждый вычислительный узел, в свою очередь, оборудован четырьмя 4-ядерными процессорами AMD Opteron 8356 2.3GHz и 32-мя гигабайтами оперативной памяти. Используемые в статье оценки производительности сетевой инфраструктуры вычислительной системы ЮФУ получены в результате выполнения тестирования в пакете Pallas MPI benchmark.

Поля скоростей водного потока, рассчитанные на основе математической модели (12)-(14), относятся к входным данным для модели транспорта вещества (1)-(2). На основе разработанных алгоритмов был построен комплекс программ, предназначенный для моделирования возможных сценариев развития экосистемы Азовского моря. На основе разработанного комплекса программ был поставлен численный эксперимент, результаты которого приведены на рис. 3-5. При численном решении задачи в начальный момент времени поле концентрации загрязняющих веществ распределено равномерно в некоторой части водоема, данная область показана на рис. 3. После чего была решена задача транспорта веществ, результаты расчетов приведены на рис. 4. Палитрой показано значение концентрации вещества.

Л - -

rrrJ-^^' -

- и - - - - г : ; - ? и

/tir i i - ;; ib: 1 i i 1 - loZ^ Ч-!-

ШштЩ:

: : : Г ? = § Н ШЩ ' ущ

® 1 = пИпПП?!? 1J

% Jii^HiiU

Рис. 3. Исходное распределение концентрации загрязняющих веществ (начальное распределение для прямой задачи)

Максимальное значение концентрации загрязняющих веществ для решенной прямой задачи (рис. 5) составляет 45 % от исходного значения концентрации в точке выброса (рис. 4). Максимальное значение концентрации загрязняющих веществ для обратной задачи (рис. 6) составляет 65 %.

Рис. 4. Распределение концентрации загрязняющих веществ по окончании заданного интервала времени (численное решение прямой задачи)

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

Рис. 5. Восстановленное поле концентрации загрязняющих веществ (решение обратной задачи)

Ниже приведены результаты численных экспериментов по моделированию течений в Азовском море (рис. 6). Палитрой показана интенсивность течения.

Из рис. 3-5 видно, что при решении обратной задачи транспорта вещества (рис. 5) исходное поле концентрации (рис. 3) восстанавливается частично, однако, несмотря на это, предложенная математическая модель позволяет достаточно точно восстановить сценарий распространения загрязняющих веществ (что является обобщением задачи проверки сценария распространения загрязняющих веществ [15]). Из приве-

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

Рис. 6. Поле скорости движения водной среды при западном ветре 5 м/с (баротропные течения)

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

БИБЛИОГРАФИЧЕСКИЙ СПИСОК

1. Самарский А.А., Вабищевич П.Н., Васильев В.И. Итерационное решение ретроспективной обратной задачи теплопроводности // Математическое моделирование. - 1997. - №. 9:5. - C. 119-127.

2. ТихоновА.Н., Арсенин В.Я. Методы решения некорректных задач. - М.: Наука, 1979.

3. Самарский А.А., Вабищевич П.Н. Численные методы решения обратных задач математической физики. - М.: Эдиториал УРСС, 2004. - 784 ^

4. СухиновА.И., ЧистяковА.Е., ТимофееваЕ.Ф., ШишеняА.В. Математическая модель расчета прибрежных волновых процессов // Математическое моделирование. - 2012. - Т. 24, № 8. - С. 32-44.

142 c.

5. СухиновА.И., НикитинаА.В., ЧистяковА.Е. Моделирование сценария биологической реабилитации Азовского моря // Математическое моделирование. - 2012. - Т. 24, № 9.

- С. 3-21.

6. Сухинов А.И., Чистяков А.Е. Параллельная реализация трехмерной модели гидродинамики мелководных водоемов на супервычислительной системе // Вычислительные методы и программирование: Новые вычислительные технологии. - 2012. - Т. 13. - С. 290-297.

7. Сухинов А.И., Чистяков А.Е., Алексеенко Е.В. Численная реализация трехмерной модели гидродинамики для мелководных водоемов на супервычислительной системе // Математическое моделирование. - 2011. - Т. 23, № 3. - С. 32-21.

8. Белоцерковский О.М., Гущин В.А., Щенников В.В. Метод расщепления в применении к решению задач динамики вязкой несжимаемой жидкости // Вычислительная математика и математическая физика. - 1975. - № 15:1. - С. 197-207.

9. Белоцерковский О.М. Турбулентность: новые подходы. - М.: Наука, 2003. - 286 с.

10. Самарский А.А. Теория разностных схем. - М.: Наука, 1989. - 616 с.

11. Самарский А.А., Николаев Е.С. Методы решения сеточных уравнений. - М.: Наука, 1978. - 592 с.

12. КоноваловА.Н. К теории попеременно-треугольного итерационного метода // Сибирский математический журнал. - 2002. - № 43:3. - С. 552-572.

13. Сухинов А.И., Чистяков А.Е. Адаптивный модифицированный попеременно-треугольный итерационный метод для решения сеточных уравнений с несамосопряженным оператором // Математическое моделирование. - 2012. - Т. 24, № 1. - С. 3-21.

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

- С. 134-154.

15. СухиновА.А. Реконструкция экологической катастрофы в Азовском море на основе математических моделей // Математическое моделирование. - 2008. - Т. 20, № 6. - С. 15-22.

REFERENCES

1. Samarskiy A.A., Vabishchevich P.N., Vasil'ev V.I. Iteratsionnoe reshenie retrospektivnoy obratnoy zadachi teploprovodnosti [Iterative solution of a retrospective inverse heat conduction problem], Matematicheskoe modelirovanie [Mathematical modeling], 1997, No. 9:5, pp. 119-127.

2. Tikhonov A.N., Arsenin V.Ya. Metody resheniya nekorrektnykh zadach [Methods of solving ill-posed problems]. Moscow: Nauka, 1979, 142 p.

3. Samarskiy A.A., Vabishchevich P.N. Chislennye metody resheniya obratnykh zadach matematicheskoy fiziki [Numerical methods for solving inverse problems of mathematical physics]. Moscow: Editorial URSS, 2004, 784 p.

4. Sukhinov A.I., Chistyakov A.E., Timofeeva E.F., Shishenya A.V. Matematicheskaya model' rascheta pribrezhnykh volnovykh protsessov [The mathematical model of calculation of coastal wave processes], Matematicheskoe modelirovanie [Mathematical modeling], 2012, Vol. 24, No. 8, pp. 32-44.

5. Sukhinov A.I., Nikitina A.V., Chistyakov A.E. Modelirovanie stsenariya biologicheskoy reabilitatsii Azovskogo morya [The scenario modeling biological rehabilitation of the Azov sea], Matematicheskoe modelirovanie [Mathematical modeling], 2012, Vol. 24, No. 9, pp. 3-21.

6. Sukhinov A.I., Chistyakov A.E. Parallel'naya realizatsiya trekhmernoy modeli gidrodina-miki melkovodnykh vodoemov na supervychislitel'noy sisteme [Parallel implementation of the three-dimensional model of hydrodynamics of shallow reservoirs on supercomputing system], Vychislitel'nye metody i programmirovanie: Novye vychislitel'nye tekhnologii [Computational methods and programming: New computational technologies], 2012, Vol. 13, pp. 290-297.

7. Sukhinov A.I., Chistyakov A.E., Alekseenko E.V. Chislennaya realizatsiya trekhmernoy modeli gidrodinamiki dlya melkovodnykh vodoemov na supervychislitel'noy sisteme [Numerical realization of a three-dimensional model of hydrodynamics for shallow reservoirs on supercomputing system], Matematicheskoe modelirovanie [Mathematical modeling], 2011, Vol. 23, No. 3, pp. 32-21.

8. Belotserkovskiy O.M., Gushchin V.A., Shchennikov V.V. Metod rasshchepleniya v primenenii k resheniyu zadach dinamiki vyazkoy neszhimaemoy zhidkosti [The method of splitting the application to solving problems of the dynamics of a viscous incompressible fluid], Vychislitel'naya matematika i matematicheskaya fizika [Computational mathematics and mathematical physics], 1975, No. 15:1, pp. 197-207.

9. Belotserkovskiy O.M. Turbulentnost': novye podkhody [Turbulence: new approaches]. Moscow: Nauka, 2003, 286 p.

10. Belotserkovskiy O.M. Turbulentnost': novye podkhody [The theory of difference schemes]. Moscow: Nauka, 2003, 286 p.

11. Samarskiy A.A., Nikolaev E.S. Metody resheniya setochnykh uravneniy [Methods for solving grid equations]. Moscow: Nauka, 1978, 592 p.

12. Konovalov A.N. K teorii poperemenno-treugol'nogo iteratsionnogo metoda [The theory of alternating triangular iterative method], Sibirskiy matematicheskiy zhurnal [Siberian mathematical journal], 2002, No. 43:3, pp. 552-572.

13. Sukhinov A.I., Chistyakov A.E. Adaptivnyy modifitsirovannyy poperemenno-treugol'nyy iteratsionnyy metod dlya resheniya setochnykh uravneniy s nesamosopryazhennym operatorom [Adaptive modified alternating triangular iterative method for solving grid equations with non-self-adjoint operator], Matematicheskoe modelirovanie [Mathematical modeling], 2012, Vol. 24, No. 1, pp. 3-21.

14. Voevodin V.V., Voevodin Vl.V. Parallel'nye vychisleniya [Parallel computing]. St. Petersburg: BKhV-Peterburg, 2004, pp. 134-154.

15. Sukhinov A.A. Rekonstruktsiya ekologicheskoy katastrofy v Azovskom more na osnove matematicheskikh modeley // Matematicheskoe modelirovanie [Mathematical modeling], 2008, Vol. 20, No. 6, pp. 15-22.

Статью рекомендовал к опубликованию д.ф.-м.н., профессор А.А. Илюхин.

Лапин Дмитрий Вадимович - Южный федеральный университет; e-mail: dmitri.lapin@gmail.com; 347928, г. Таганрог, пер. Некрасовский, 44; тел.: 88634371606; кафедра математического обеспечения суперкомпьютеров.

Чистяков Александр Евгеньевич - e-mail: cheese_05@mail.ru; кафедра математического обеспечения суперкомпьютеров; доцент.

Сухинов Антон Александрович - e-mail: soukhinov@gmail.com. 347932, г. Таганрог, ул. Пархоменко, 58/1; тел.: 89296796206; научный сотрудник.

Lapin Dmitry Vadimovich - Southern Federal University; e-mail: dmitri.lapin@gmail.com; 44, Nekrasovskiy, Taganrog, 347928, Russia; phone: +78634371606; the department of mathematical software supercomputers.

Chistyakov Alexander Evgenjevich - e-mail: cheese_05@mail.ru; the department of mathematical software supercomputers; associate professor.

Sukhinov Anton Alexandrovich - e-mail: soukhinov@gmail.com; 58/1, Parkhomenko, Taganrog, 347932, Russia; phone: +79296796206; research scientist.

УДК 519.684.6

М.Д. Чекина

МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ЗАДАЧ ГЕОФИЛЬТРАЦИИ В ПОЧВОГРУНТАХ С ФРАКТАЛЬНОЙ СТРУКТУРОЙ НА МНОГОПРОЦЕССОРНЫХ ВЫЧИСЛИТЕЛЬНЫХ СИСТЕМАХ*

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

модели на основе уравнения Букингема-Ричардса, в котором оператор дифференцирования по

*

Работа выполнена при поддержке Министерства образования и науки РФ по Соглашению с уникальным идентификатором РГМЕП57814Х0006.