автореферат диссертации по информатике, вычислительной технике и управлению, 05.13.18, диссертация на тему:Математическое моделирование процессов переноса примесей в жидкостях и пористых средах

На правах рукописи

Сухинов Антон Александрович

МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ПРОЦЕССОВ ПЕРЕНОСА ПРИМЕСЕЙ В ЖИДКОСТЯХ И ПОРИСТЫХ СРЕДАХ

Специальность 05.13.18 —Математическое моделирование, численные методы и комплексы программ

Автореферат ? 2 л0.п ?гг{?

диссертации на соискание учёной степени кандидата физико-математических наук

Москва - 2009

003482909

Работа выполнена в Институте математического моделирования Российской академии наук

Научный руководитель: доктор физико-математических наук, член-

корреспондент РАН, профессор Четверушкин Борис Николаевич

Официальные оппоненты: доктор физико-математических наук, профес-

сор Мажукин Владимир Иванович

доктор физико-математических наук, профессор Шевелёв Юрий Дмитриевич

Ведущая организация: Институт вычислительной математики РАН,

г. Москва

Защита состоится « 3 » декабря 2009 года в 14 ч. 20 мин. на заседании Диссертационного совета Д.002.058.01 при Институте математического моделирования РАН по адресу: 125047, Москва, Миусская площадь, 4а.

С диссертацией можно ознакомиться в библиотеке Института математического моделирования РАН.

Автореферат разослан * » октября 2009 года.

Учёный секретарь

Диссертационного совета Д.002.058.01 доктор физико-математических наук

I/

Н.В. Змитренко

1 Общая характеристика работы

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

1.1 Актуальность работы

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

Численное решение задач математической физики обычно требует наличия сетки, покрывающей расчётную область. Дискретизация уравнений модели производится на уровне элементов сетки. Физические явления часто характеризуются наличием разрывов и локализованных областей с высокими градиентами физических величин. Численное моделирование таких задач и используемые методы (конечные элементы, конечные объёмы, конечные разности и спектральные методы) на равномерных сетках неэффективны, когда требуется высокая точность решения. Точность может быть повышена путём использования схем высокого порядка аппроксимации или путём уменьшения размеров элементов сетки. Первый метод имеет недостаток: схемы высокого порядка аппроксимации приводят к появлению осцилляций [1]. Уменьшение же размеров ячеек во всей области значительно увеличивает вычислительную сложность задачи.

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

объёма диссертации посвящена адаптивным сеткам).

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

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

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

1.2 Цели и задачи работы

Целями работы являются:

• Разработка эффективных алгоритмов и компьютерных программ построения декартовых сеток, динамически адаптирующихся к решению.

• Применение разработанных методов к решению задачи переноса примесей в пористой среде.

• Разработка и исследование модели переноса многокомпонентной примеси в водной среде.

• Применение модели переноса многокомпонентной примеси для моделирования процесса возникновения анаэробной зоны в Азовском море.

• Разработка алгоритмов обработки исходных данных: измерение течений, концентрации примесей, формы донной поверхности.

1.3 Научная новизна

Научная новизна работы заключается в новых методах и моделях:

• Комплекс алгоритмов адаптации декартовых сеток, имеющий низкие накладные

расходы.

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

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

• Метод интерполяции пространственно распределённых двумерных данных, учитывающий форму области интерполяции, погрешность исходных данных, и наличие поля выделенных направлений (например, течений в водоёме).

Метод адаптации декартовых сеток имеет следующие отличия от аналогичных методов:

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

• Результат адаптации сетки не зависит от порядка обработки ячеек. Это позволяет избежать «зацикливания» алгоритма адаптации и минимизировать изменение сеточной функции, вносимое адаптацией.

• Алгоритм адаптации и генерируемые им структуры сеток устойчивы относительно сеточной функции и их изменения во времени.

1.4 Теоретическая и практическая значимость

Разработанные алгоритмы адаптивных декартовых сеток могут быть использованы при численном решении уравнений математической физики различной размерности. Некоторые из разработанных методов (например, метод интерполяции сеточной функции и построения разностной схемы, методы решения СЛАУ на адаптивных сетках) могут быть применены независимо от остального комплекса алгоритмов.

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

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

Все рассмотренные в диссертационной работе алгоритмы реализованы в виде компьютерных программ на языках программирования Delphi и С++.

1.5 На защиту выносятся:

• Алгоритмы декартовых сеток с динамической адаптацией к решению, в том числе:

- структуры данных, необходимые для хранения сеток и их эффективной адаптации к получаемому решению;

- алгоритм адаптации сетки;

- методы интерполяции сеточной функции и аппроксимации дифференциальных уравнений на таких сетках;

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

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

• Метод интерполяции исходных данных, учитывающий погрешность этих данных, наличие поля выделенных направлений, и форму области интерполяции

1.6 Апробация работы:

Изложенные в диссертации результаты докладывались на 10-й Международной конференции «Mathematical Modeling and Analysis» (Тракай, Литва, 2005), на конференции «European Conference on Computational Fluid Dynamics» (Егмон Аан Зи, Нидерланды, 2006), третьей Международной конференции «Adaptive Modeling and Simulation» (Гётеборг, Швеция, 2007), Международной конференции «Parallel Computational Fluid Dynamics» (Анталия, Турция, 2007), Научно-практическом семинаре «Круглый стол МФТИ-Шлюмберже» (Москва, 2007), 50-й Научной конференции МФТИ «Современные проблемы фундаментальных и прикладных наук» (Долгопрудный, 2007), третьей Международной конференции «Computational Methods in Applied Mathematics» (Минск, Белоруссия, 2007), пятой Всероссийской межву-

зовской конференции молодых учёных (Санкт-Петербург, 2008), конференции «Параллельные вычислительные технологии» (Нижний Новгород, 2009).

1.7 Структура и объём работы:

Диссертация состоит из введения, четырёх глав, заключения и списка литературы. Работа изложена на 150 страницах, включает 35 рисунков и 5 таблиц. Список литературы включает 40 наименований.

1.8 Публикации по теме диссертации:

По теме диссертационной работы опубликовано 7 статей. Их список представлен в конце автореферата.

2 Содержание работы

2.1 Введение

Введение включает обоснование актуальности темы, содержит формулировку основных целей работы и краткое содержание глав.

2.2 Первая глава

В первой главе приводится классификация существующих алгоритмов адаптивных сеток, и описываются отличия предлагаемых методов от оригинальной работы "Adaptive mesh refinement for hyperbolic partial differential equations" [2].

Затем вводятся основные предположения адаптивных сеток (на примере двумерных сеток, хотя алгоритмы легко обобщаются на пространство любой размерности):

1) Корневая ячейка сетки: вначале сетка состоит из единственной прямоугольной ячейки R С R2, которая называется корневой ячейкой сетки.

2) Разбиение ячеек: ячейка может быть разбита (разделена) на четыре ячейки одинакового размера. Эти меньшие ячейки называются дочерними ячейками (потомками) разделённой ячейки. Сама разделяемая ячейка называется родителем своих потомков.

3) Объединение ячеек: четыре ячейки могут быть объединены в одну, только если они когда-то составляли одну ячейку.

4) Сеточная функция: каждой ячейке сетки С поставлена в соответствие величина / (С) € К" (скалярная или векторная), описывающая среднее значение вычис-

ляемой функции в пределах ячейки. Величину / (С) будем называть значением ячейки С.

5) Ограничение на размеры ячеек: ячейка не должна иметь в своей окрестности ячеек, которые отличаются от неё по размерам более чем в два раза.

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

Рисунок 1: Четверичное дерево — удобный способ представления двумерной адаптивной сетки

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

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

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

Определяется процедура сшивки адаптивных сеток, а также методы построения фиктивных сеток для периодических граничных условий, граничных условий на непроницаемой границе, и граничных условий первого рода.

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

Ключевое место в алгоритмах адаптивных сеток занимают алгоритмы интерполяции сеточной функции. Интерполяция используется как для построения разностных схем (например, методом конечных объёмов), так и для оценки «качества» ячеек при адаптации сетки.

Имея значения ячеек, можно построить непрерывно дифференцируемую функцию /(х,у), которая удовлетворяет следующим уравнениям для каждой ячейки С:

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

К сожалению, построение такой функции обычно требует решения системы линейных уравнений. Поэтому вместо непрерывно дифференцируемой функции / построим непрерывные (но не непрерывно дифференцируемые) функции /, /х и /у, которые аппроксимируют /, и /у соответственно. При аппроксимации вторых производных методом конечных объёмов будем использовать /х и /у вместо дифференцирования функции /.

Разобьём ячейку на 4 интерполяционных прямоугольника Р\, Р-х, Рз, Р4, построенных по девяти интерполяционным точкам ро,...,ря■ Каждая точка Pi (г = = 0,..., 8) содержит значения /,, которые аппроксимируют неизвестные зна-

(1)

с

Рк

Р7

Р5

—у 1 1

Р1 1 Р2 1

_______, 1

Рз ! Р4 1

• —А......

РЗ

Р8

Рб

Р4

Рисунок 2: Малые соседи и интерполяционные точки

чения /х (р{) и /у(рг) соответственно. В пределах каждого прямоугольни-

ка величины йй- интерполируются билинейным способом, формируя функции /(х, у), /х (х, у) и р (х, у).

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

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

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

мера друг относительно друга.

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

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

листовой ячейки С* определяется следующим выражением:

Д = (2)

В скалярном случае (когда /? — скалярные величины), вариация —это разность между максимальным и минимальным значениями интерполяционной функции в пределах ячейки:

Д = max (/Л - min (/j). (3)

Цель алгоритма адаптации сетки — минимизация максимальной вариации среди листовых ячеек сетки при ограниченном числе этих ячеек:

max (Di) -> min при iVleafs < Nneftä, (4)

i 61

где L — множество листовых ячеек, iViea[S — число листовых ячеек, iVneed ~ максимально разрешённое число листовых ячеек в сетке.

Интерполяционная функция / фиксируется перед изменением топологии сетки (не перестраивается в процессе адаптации). Это важно для следующего:

1. Нет необходимости пересчитывать интерполяционную функцию после каждого разбиения/объединения во время процесса перестроения.

2. Значения ячеек (вычисляемые через интерполяционную функцию при разбиении ячеек) не будут зависеть от способа, которым текущее состояние сетки получено

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

В работе приведён и исследован алгоритм адаптации сетки, основная идея которого—разбиение ячеек с максимальной вариацией и объединение ячеек с минимальной вариацией. Доказано, что алгоритм обладает следующими свойствами, если сеточная функция не равна константе во всей сетке (иначе вся сетка может быть представлена единственной ячейкой), и если не установлены запреты на разбиение или объединение некоторых ячеек (однако могут быть установлены ограничения на максимальный и минимальный размер ячеек):

1. Алгоритм перестроения минимизирует максимальную вариацию данных в сетке (он достигает цели (4)).

2. Если не требуется изменение числа ячеек в сетке (./V = Агпее(1 « Л^а^), то алгоритм адаптации работает за время 0(Щ в среднем, и О (ЛЧп(Лг)) в худшем случае.

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

Для решения этих проблем в работе предложены специальные методы упорядочивания узлов в соответствии с кривой Гильберта (рисунок 3, сверху). Такой обход узлов позволяет увеличить скорость сходимости и расширить область сходимости итерационных методов по сравнению с широко используемыми способами обхода дерева сетки в ширину и в глубину. Кроме того, предложен метод «раскраски» ячеек адаптивной сетки в 8 цветов, являющийся обобщением красно-чёрного упорядочивания узлов на равномерной сетке.

Оригинальной идеей, изложенной в работе, является смена порядка обхода узлов после каждой итерации решения СЛАУ (рисунок 3, снизу). Такое изменение порядка обхода позволяет добиться дальнейшего улучшения сходимости. На рисун-

Fif T--F1

-МЧ F--r T--F-

Ш Wf \ш

щ №fl

±\ ftt

г £ Г1

1 ч- 1

Рисунок 3: Сверху — построение кривой Гильберта при измельчении сетки. Снизу — поворот кривой после каждой итерации решения СЛАУ

ке 4 показана зависимость требуемого количества итераций от числа Куранта при решении модельной задачи конвекции на сетке, состоящей из 4096 ячеек.

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

Первая глава заканчивается проверкой разработанных методов на простейшей двумерной задаче конвективного переноса:

дс дс дс

(5)

где с = с(х, у, t) — некоторая физическая величина, £ — время, (vx,vy) = const — вектор скорости переноса величины с.

Область решения — единичный квадрат, граничные условия — с = 0 на границе области, начальное значение — с (х, у, 0) = со (х, у).

Несмотря на простоту постановки задачи и известное аналитическое решение, она сложна для численной реализации сеточными методами. Если при разностной аппроксимации задачи не внести искусственную диффузию («схемную вязкость»), то разностная схема будет неустойчива. Если же схемную вязкость внести (например, путём использования разностной схемы, ориентированной против потока), то численное решение будет иметь большую погрешность, так как в исходном диф-

Рисунок 4: Зависимость требуемого числа итераций метода Гаусса-Зейделя от числа Куранта для различных обходов элементов сетки

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

Пусть начальные условия —это пятно единичной концентрации в форме круга диаметром 0,25, центр которого имеет координаты (0,25;0,25). Поле скоростей — вектор (1; 1) во всей области. Через время 0,5 центр пятна должен находиться в точке (0,75; 0,75).

На рисунке 5 показаны результаты вычислений на равномерной и адаптивной сетках. Обе сетки имеют по 4096 ячеек. Можно видеть, что результат, полученный на адаптивной сетке, гораздо меньше «размазан». Круг имеет внутри себя большую область с почти постоянным значением искомой функции. В то же время результат, полученный на обычной сетке, размазан значительно. Точность полученных решений приведена в таблице 1.

Рисунок 5: Результаты вычисления тестовой задачи на равномерной (слева) и адаптивной (справа) сетках. Сетки состоят из 4096 ячеек каждая. Пунктиром показано начальное положение пятна

Таблица 1: Среднеквадратичное отклонение результатов моделирования

от точного решения

2.3 Вторая глава

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

Нефтяной пласт предполагается тонким, поэтому задача рассматривается в двумерной постановке. Схемы расположения скважин могут быть различными. Характерные расстояния в задаче — сотни и тысячи метров, характерные времена — месяцы и годы. В этих предположениях отношение капиллярного давления к полному

1024 ячейки 4096 ячеек 16384 ячейки

Равномерная сетка 0,125 Адаптивная сетка 0,07

0,09 0,04

0,07 0,035

гидродинамическому падению давления мало. Это позволяет пренебречь капиллярными силами, полагая, что течение двух фаз: т (воды) и о (нефти) в пористой среде подчиняется классической модели Бакли-Леверетта. Модель основана на следующих предположениях: жидкости несмешиваемы и несжимаемы, среда недефор-мируема, течения фаз подчиняются закону Дарси, капиллярные и гравитационные силы не учитываются. Пористость в области мы считаем постоянной; проницаемость — неоднородной.

Модель Бакли-Леверетта двухфазной фильтрации в пористой среде описывается следующими уравнениями:

яе I а ■ Р(5"*) — на источниках;

+ <Иу (Р (5) К (5) ега<1Р) = <Г = < (6)

I <7 • Р (й1) — в остальной области;

<11у(К(5)ёгас1Р) =9. (7)

Здесь т — пористость; 5 — водонасыщенность; Р (.5?) — функция Бакли-Леверетта; К (5) — нелинейный коэффициент, включающий проницаемости и вязкости; Р — давление; q — дебиты на скважинах; 5* — критическая водонасыщенность.

Скорость фильтрации флюида (закон Дарси): = К (й) gradP; скорость фильтрации воды: \¥,ц = Р(5)\У/. Нелинейные коэффициенты:

Здесь кш (5) — относительная проницаемость воды; к0 (¿>) — относительная проницаемость нефти; к — абсолютная проницаемость; цш - вязкость воды; р0 — вязкость нефти.

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

+ а) + (Цу + Э О*,)) = Яс = • С*. (9)

Здесь С — концентрация примеси в воде, а — концентрация примеси, адсорбированной на стенках пор (а = а (С) в случае обратимой адсорбции), фс — источник

Рисунок 6: Рассчитанные поля водонасьнценности и концентрации примеси. Каждая сетка имеет 4096 ячеек

примеси, С* — концентрация примеси в источнике, S — поток, вызванный конвективной диффузией:

дС

Si = -Dij—, i,je{x,y}. (10)

Dij — эффективный тензор конвективной диффузии. Имеются феноменологические формулы для его получения из скорости Дарси и некоторых характеристик среды. В расчётах использована формула В.Н. Николаевского [4]:

Аз = [(Ai - А2) 5ij + \2щпj) Ww, (11)

где Ww = |W„,|, n = Ww/Ww, Xi и Лг — некоторые положительные коэффициенты. Можно видеть, что даже в изотропной среде имеется выделенное направление, определяемое вектором скорости Дарси.

Проведены расчёты для пятиточечной тестовой задачи с использованием адаптивных сеток. Физическая область квадратная (150 х 150 метров). Нагнетающая скважина располагается в юго-западном углу области, продуктивная скважина — в северо-восточном. Поле проницаемости существенно неоднородное: абсолютная проницаемость меняется в пределах от 10~8 до Ю-12 м2. Другие параметры модели: т = 0,2, дебит скважины: 300 м3/сутки, толщина нефтяного пласта: 10 метров (нужна для перевода дебита в двумерную постановку), цт = 1,2 • 10"13, ц0 =

= 6- Ю-13. Результаты моделирования (100 дней добычи) показаны на рисунке 6. Максимумы водонасьпценности, концентрации и давления соответствуют водона-гнетающей скважине, через которую в пласт поступает примесь. Движущиеся фронты имеют сложную форму (образуются «языки»). Для поля давления было выделено в 4 раза меньше ячеек, чем для других величин. Это позволило значительно снизить время вычислений, так как задача расчёта давления наиболее трудоёмкая.

2.4 Третья глава

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

^ + + = + г = (12)

где С1 — концентрация примеси г; ¿ — время; V = I д/дх + j д/ду + к д/дг — оператор набла; V : (У,у) = 0 —скорость течения жидкости; gi —скорость гравитационного осаждения примеси г; Ф : (х, Фх) = 0 — тензор диффузии примеси г; N

Т : У; = 1,..., Ы, ^ = 0 — матрица взаимного преобразования примесей; — 1=1

источники/стоки примеси г.

Модель применяется для исследования возможной причины возникновения анаэробной зоны, которая была обнаружена в августе 2001 года во время экспедиции в Азовском море [5].

Зона аноксии располагалась в восточной части моря (рисунок 7). В этом районе наблюдалась устойчивая стратификация водной среды. Придонный слой воды толщиной от 1 до 4-х метров характеризовался увеличенной плотностью, отсутствием кислорода и высоким содержанием сероводорода. Он отделялся от верхних слоёв переходной зоной (толщиной 0,5-1 м), максимальный градиент плотности в которой составлял 10% на метр. На границе придонного слоя также наблюдался слой мутности, характерный для водоёмов с анаэробными условиями (Чёрное море, Балтийское море). Концентрация сероводорода на глубине 7-8 метров превышала его концентрацию в Чёрном море на глубинах более 350 метров. В результате этого явления произошел массовый замор рыбы.

Автор воспользовался моделью течений, разработанной для Азовского моря другими исследователями [6], чтобы проверить гипотезу о природном возникновении

I Рисунок 7: Концентрация растворённого кислорода (мкм) в придонном слое. [ Пунктирная линия показывает границу области с анаэробными условиями

I

' этого явления [7]. Гипотеза заключается в том, что высокое содержание сероводорода (и низкое содержание кислорода) в определённом участке моря вызвано большим ' количеством органических примесей, попавших в воду вместе со стоками рек. Эти

примеси выделяли сероводород в процессе своего разложения. I Модель переноса примесей была адаптирована к требуемой задаче. Результат мо-I делирования показан на рисунке 8.

2.5 Четвёртая глава

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

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

' Основное внимание в этой главе уделено методу интерполяции исходных данных, разработанному и исследованному автором. Метод является разновидностью ! метода наименьшей кривизны [8]. Математически задача формулируется, как оп-

Рисунок 8: Рассчитанная зона аноксии. Числа показывают концентрацию кислорода (мкм)

тимизационная. Для получения результата должна быть решена система линейных уравнений с ограничениями. Это требует большого объёма вычислений при использовании итерационных методов. Для существенного сокращения количества вычислений использован многомасштабный подход.

Применительно к задаче интерполяции глубины метод формулируется следующим образом. В плоской области G заданы ограничения на диапазон изменения глубин: функции т(х,у) и М(х,у) такие, что т(х,у) < М(х,у), где (х,у) G 6 G — горизонтальные координаты некоторой точки области. Это означает, что значение реальной глубины/высоты е(х,у) в области G удовлетворяет неравенствам т(х,у) ^ е(х,у) ^ М(х,у). Допустимы случаи, когда в некоторых точках или подобластях тп(х,у) = -оо и/или М (х,у) = +оо.

Требуется построить в области G функцию f (х,у), которая имеет непрерывные производные первого порядка и ограниченные производные второго порядка, удо-

влетворяющую следующим условиям:

т{х,у)^!{х,у)^М(х,у), (13)

//Н2МЗ шнвишу-^^

(х,з/)€С 4 '

где а > 0 — параметр, значение которого будет описано ниже. Постановка задачи имеет следующие отличия от классической задачи интерполяции [8]:

• исходные данные представлены в интервальном виде: [т (х, у), М (х, у)];

• форма области <3 также берётся в расчёт; это важно при интерполяции физических величин, существующих только в пределах некоторой области (например, концентрации загрязнений в водоёме).

Перечисленные особенности расширяют область применения описываемого интерполяционного метода.

Минимизация суммы квадратов первых производных даёт тенденцию горизонтального положения интерполированной поверхности. Это соответствует минимизации потенциальной энергии земной поверхности. Коэффициент а ^ 0 (обычно достаточно малый) предназначен для настройки этой особенности; его оптимальное значение может быть оценено на основе анализа исходных данных.

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

В качестве примера рассмотрена задача интерполяции данных, показанных на рисунке 9. Следует отметить, что метод не требует распознавания изолиний глубины: в качестве входных данных могут быть использованы цвета областей, поэтому достаточно залить области между изолиниями различными цветами. Каждый цвет обозначает свой диапазон глубин, который записывается в ограничения по высоте при построении сетки. Кроме областей с цветом карта имела множество отметок глубины и высоты (на рисунке 9 показаны кружочками), которые были оцифрованы и внесены в соответствующие ячейки сетки при её построении. Погрешность всех данных положена равной 0,1 метра. Кроме того, была отдельная карта для треугольного водоёма на юго-востоке карты (Этанг-де-Больмон). Её данные не показаны.

Для вычислений использовалась сетка размерами 512 х 512 ячеек. Вначале был произведён расчёт одномасштабным алгоритмом. Потребовалось более 50-ти тысяч итераций (14 минут машинного времени) для того, чтобы итерации сошлись.

5 км.

Данные о высоте

1 о < н

| -3 м. < Н < О -5 м. 5 Н < -3 м. Н < -5 м.

Высота неизвестна О Точки с заданной высотой

Рисунок 9: Пример данных о высоте для залива Этанг-де-Берр (Франция)

Затем был произведён расчёт многомасштабным алгоритмом. На сетке каждого масштаба требовалось произвести от 10 до 100 итераций. В итоге алгоритм сошёлся за 215 итераций (2,9 секунд машинного времени). Использование многомасштабного алгоритма позволило ускорить вычисления в 300 раз.

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

2.6 Основные результаты и выводы

• Разработаны эффективные алгоритмы динамических адаптивных сеток.

5000

10000

15000

20000

Рисунок 10: Результат интерполяции данных с рисунка 9.

Сетка 512 х 512 ячеек, а = 3 - 10~3 1/м. Масштаб указан в метрах

• Алгоритмы реализованы в виде библиотеки компьютерных программ и опробованы на задаче переноса пассивной примеси в пористой среде

• Разработана и исследована модель переноса многокомпонентной примеси в водной среде.

• На основе разработанной модели подтверждена гипотеза о механизме возникновения анаэробной зоны в Азовском море.

• Разработаны алгоритмы обработки и интерполяции исходных данных.

2.7 Основные результаты опубликованы в работах:

1. A.I. Sukhinov, В. Roux, A.A. Sukhinov Reconstruction of Basin Bottom Surface for High Precision Hydrodynamics Modeling Using Parallel Computations. Int. Conf on Parallel Computational Fluid Dynamics, May 21-24, 2007, Antalya, Turkey (Parallel CFD 2007). CD-ROM Proceedings, ParCFD-2007-031, 8 pp.

2. B.N. Chetverushkin, N.G. Churbanova, M.A. Trapeznikova, A.A. Sukhinov, A.A. Malinovskij, Adaptive Cartesian mesh refinement for simulating multiphase flows in porous media.

Computational Methods in Applied Mathematics, Vol. 8 (2008), No 2, Minsk, Belarus, pp. 101-115.

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

4. А.А. Сухинов, Адаптивные декартовы сетки. Научно-технический вестник Санкт-Петербургского государственного университета информационных технологий, механики и оптики, № 54 (2008), стр. 164-165.

5. А.А. Сухинов, Восстановление донной поверхности по различным картографическим данным. Вестник Южно-Уральского Государственного Университета, вып. 3 (2009), 10 стр.

6. Boris N. Chetverushkin, Natalia G. Churbanova, Anton A. Sukhinov, Marina A. Trapezni-kova, Technique of Cartesian Mesh Refinement with Dynamic Adaptation to the Solution. Adaptive Modeling and Simulation 2009. Proc. oí the IV International Conference on Adaptive Modeling and Simulation, Bruxelles, Belgium, 25-27 May, 2009, ed. by Ph. Bouillard and P. Diez, CIMNE, Barcelona, pp. 55-58.

7. A.A. Sukhinov, Mathematical Reconstruction of the Ecological Disaster in the Sea of Azov. Mathematical Models and Computer Simulations, 2009, Vol. 1, No. 3, pp. 353-359.

2.8 Цитируемая литература

1. С.К. Годунов, Разностный метод численного расчёта разрывных решений уравнений гидродинамики. Мат. сборник, 1959, т. 47(89):3. Стр. 271-306.

2. M.J. Berger, J. Oliger, Adaptive mesh refinement for hyperbolic partial differential equations. Journal of Computational Physics. Vol. 53, pp. 484-512. Mar. 1984.

3. R.A. Finkel, J.L. Bentley, Quad trees: A Data Structure for Retrieval on Composite Keys. Acta Informática 4(1), pp. 1-9, 1974.

4. B.H. Николаевский, Механика пористых и трещиноватых сред, Москва, Недра, 1984.

5. Е.В. Якушев, А.И. Сухинов и др. (2003). Комплексные океанологические исследования в Азовском море во время 28-й экспедиции научно-исследовательского судна "Акванавт" (Июль-Август 2001 г.). Журнал «Океанология», РАН, т. 43, стр. 44-53.

6. В.С.Васильев, А.И.Сухинов (2003). Прецизионные двумерные модели мелководных водоёмов. «Математическое моделирование», РАН, т. 15, № 10, стр. 17-34.

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

8. Schu, G.H. Review of Interpolation Methods for Digital Terrain Models. Canadian Surveyor, December, 1976. Vol. 30, N. 5. P. 389-412.

Подписано в печать 22.10.2009 г. Заказ 135. Тираж 100 экз.

Отпечатано в ООО "ФЭД+", Москва, ул. Кедрова, д. 15, тел. 774-26-96 Формат 60x90/18. Печать офсетная. Объем 1,25 п. л.

Введение

Глава 1. Ортогональные сетки с динамической адаптацией к решению

1.1. О бзор оригинальной концепции Adaptive Mesh Refinement (M.J. Berger, J. Oliger, P. Colella).

1.2. Концепция Building Cube.

1.3. Предлагаемая концепция

1.4. Основные предположения.

1.5. Структура сетки.

1.6. Сеточная функция.

1.7. Преобразование сеток.

1.8. Адаптация сетки.

1.9. Неравномерный шаг по времени.

1.10. Быстродействие алгоритмов и его повышение.

1.11. Хранение и передача данных сетки.

1.12. Итерационные методы решения СЛАУ на адаптивной сетке

1.13. Параллельный алгоритм для систем с распределённой памятью

1.14. Обобщение алгоритмов на многомерный случай.

1.15. Тестовая задача переноса.

Глава 2. Перенос пассивной примеси в пористой среде.

2.1. Основные предположения.

2.2. Уравнения фильтрации в пористой среде.

2.3. Уравнения переноса пассивной примеси.

2.4. Применение адаптивных сеток.

2.5. Результаты численных экспериментов.

Глава 3. Перенос примесей в водной среде.

3.1. Концентрация примесей и мутность.

3.2. Перенос примесей

3.3. Модель переноса многокомпонентной примеси.

3.4. Моделирование формирования анаэробной зоны в Азовском море

Глава 4. Измерение и обработка исходных данных.

4.1. Измерение мутности оптическим способом

4.2. Интерполяция исходных данных и реконструкция донной поверхности

Данная работа посвящена различным аспектам решения задач математической физики: от методов построения расчётных сеток и решения сеточных СЛАУ до обработки экспериментальных данных и верификации моделей.

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

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

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

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

Институт математического моделирования Российской академии наук

На правах рукописи

Сухинов Антон Александрович

Математическое моделирование процессов переноса примесей в жидкостях и пористых средах

05.13.18 - Математическое моделирование, численные методы и комплексы программ

ДИССЕРТАЦИЯ на соискание ученой степени кандидата физико-математических наук

Научный руководитель д. ф.-м. н., проф. Четверушкин Борис Николаевич

Москва - 2009

Содержание

Введение .4

Заключение

Разработаны новые алгоритмы двумерных адаптивных сеток. Метод адаптивных сеток значительно уменьшает время вычислений и снижает требуемый объём оперативной памяти по сравнению с равномерными сетками при той же точности получаемого решения. Разработаны алгоритмы для систем с общей и распределённой памятью. Разработаны некоторые сопутствующие алгоритмы, в том числе методы решения сеточных СЛАУ и методы хранения и передачи данных сетки.

Алгоритмы адаптивных сеток реализованы в виде библиотеки компьютерных программ на языке программирования С++ для их широкого применения в прикладных задачах.

Алгоритмы опробованы на задаче переноса пассивной примеси в пористой среде (пятиточечная тестовая задача).

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

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

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

1. С.К. Годунов. Разностный льетод численного расчёта разрывных решений уравнений гидродинамики, Мат. сборник, 1959, т. 47(89):3. Стр. 271-306.

2. M.J. Berger, J. Oliger, Adaptive mesh refinement for hyperbolic partial differential equations, Journal of Computational Physics. Vol. 53, pp. 484-512. Mar. 1984

3. M.J. Berger, P. Colella, Local Adaptive Mesh Refinement for Shock Hydrodynamics, Journal of Computational Physics, 82 (1989), pp. 65-84

4. W.F. Mitchell, A comparison of adaptive refinement techniques for elliptic problems. ACM Transactions on Mathematical Software, 15 (1989) 326-347.

5. L.F. Richardson, On the Approximate Arithmetical Solution by Finite Differences of Physical Problems Involving Differential Equations, with an Application to the Stresses in a Masonry Dam. Proc. R. Soc. Lond. A March 2, 1910

6. A.A. Самарский, Ю.Л. Попов, Разностные методы решения задач газовой динамики, Учеб. пособие: Для Вузов. 3-е изд., доп., Наука, Физматлит, М., 1992

7. P.J. Frey, Р.Ь. George, Mesh Generation: Application to Finite Elements, HERMES Science Publishing, Oxford, UK, 2000.

8. V.D. Liseikin, Grid Generation Methods, Springer-Verlag, Berlin, 1999.

9. M.J. Aftosmis, Solution adaptive Cartesian grid methods for aerodynamic flows with complex geometries, Lecture Notes for 28th CFD Lecture Series, von Karman Institute for Fluid Dynamics, Belgium, 1997, pp. 1-108.

10. J.A. Trangenstein, Multi-scale iterative techniques and adaptive mesh refinement for flow in porous media, Advances in Water Resources, 25 (2002), pp. 1175-1213.

11. P.W. Hemker, On the structure of an adaptive multi-level algorithm. BIT Numerical Mathematics, vol. 20, issue 3 (1980). pp. 289-301.

12. I.S. Kim, W.J.R. Hoefer, A local mesh refinement algorithm for the time domain-finite difference method using Maxwell's curl equations, IEEE Trans. Microwave Theory Tech., vol. 38, pp. 812-815, June 1990.

13. M. Okoniewski, E. Okoniewska, M.A. Stuchly, Three-dimensional subgridding algorithm for FDTD, IEEE Trans. Antennas Propagat., vol. 45, pp. 422-429, Mar. 1997.

14. Shun Takahashi, Takashi Ishida, Kazuhiro Nakahashi, Building-Cube Method for Incompressible Flow Simulations of Complex Geometries. Springer Berlin Heidelberg, Computational Fluid Dynamics (2008), part 25 , pp. 474-478.

15. T. Ishida, S. Takahashi, K. Nakahashi, Fast Cartesian Mesh Generation for Building-Cube Method Using Multi-Core PC, AIAA 2008-919, 46th AIAA Aerospace Sciences Meeting (2008).

16. M.J. Berger, M.J. Aftosmis, Aspects (and aspect ratios) of Cartesian mesh methods, in: Proc. of the 16th Int. Conf. on Numerical Methods in Fluid Dynamics, (6-10 July, Л998, Arcachon, France), Springer-Verlag, Heidelberg, Germany.

17. N. Hannoun, V. Alexiades, Issues in Adaptive Mesh Refinement Implementation, Sixth Mississippi State Conference on Differential Equations and Computational Simulations, Electronic Journal of Differential Equations, Conference 15 (2007), pp. 141-151.

18. R.A. Finkel, J.L. Bentley, Quad trees: A Data Structure for Retrieval on Composite Keys. Acta Informatica 4(1), pp. 1-9, 1974.

19. H.C. Кошляков, Э.Б. Глинер, M.M. Смирнов. Уравнения в частных производных математической физики. — М.: Высшая школа, 1970.

20. В.Н. Николаевский, Механика пористых и трещиноватых сред, Москва, Недра, 1984.

21. R.E. Ewing (Ed.), The mathematics of reservoir simulations, SIAM, Philadelphia, 1983.

22. B.M. Ентов, А.Ф. Зазовский, Гидродинамика процессов повышения нефтеотдачи, Москва, Недра 1989.

23. М.А. Trapeznikova, N.G. Churbanova, Simulation of multi-phase fluid filtration on parallel computers with distributed memory, Proc. of the 4th European CFD Conf. (Eds K. D. Papailiou et al.), Wiley, Chichester, 1998, Vol. 1, Pt. 2, pp. 929-934.

24. M.A. Trapeznikova, N.G. Churbanova, B.N. Chetverushkin, Parallel elliptic solvers and their application to oil recovery simulation, Grand Challenges in Computer Simulation, Proc. of Simulation Multi Conf., SCS, San Diego, CA, 2000, pp. 213-218.

25. L.H. Howell and J.B. Bell, An adaptive mesh projection method for viscous incompressible flow, SIAM J. Sci. Comput. 18 (1997), No. 4, pp. 996-1013.

26. B.N. Chetverushkin, N.G. Churbanova, A.A. Sukhinov, M.A. Trapeznikova, Simulation of contaminant transport in an oil-bearing stratum at water flooding, Proc. of ECCOMAS CFD 06 (Eds P. Wesseling, E. Onate and J. Periux), TU Delft, The Netherlands 2006.

27. E.B. Якушев, А.И. Сухинов и др, Комплексные океанологические исследования в Азовском море во время 28-й экспедиции научно-исследовательского судна "Акванавт" (Июль-Август 2001 г.). Океанология, РАН, т. 43 (2003), стр. 44-53.

28. E.I. Debolskaya, E.V. Yakushev, A.I. Sukhinov, Formation of Fish Kills and Anaerobic Conditions in the Sea of Azov. Water Resources, Vol. 32, No 2, 2005, pp. 151-162.

29. B.C. Васильев, А.И. Сухинов, Прецизионные двумерные модели мелководных водоёмов. Математическое моделирование, РАН, т. 15, № 10 (2003), стр. 17-34.

30. A.A. Сухинов, А.И. Сухинов, О единственности решения задачи диффузии-конвекции-агрегирования взвесей в водной среде. Тезисы международной конференции «Современные проблемы механики сплошной среды». Ростов-на-Дону, т. 1 (2003), стр. 180-184.

31. А.В. Колдоба, Ю.А. Повещенко, Е.А. Самарская, В.Ф. Тишкин Методы математического моделирования окружающей среды. М., Наука, 2000.

32. G.H. Schu,Review of Interpolation Methods for Digital Terrain Models. Canadian Surveyor. December, 1976. Vol. 30, N. 5. pp. 389-412.

33. M.A. Oliver, R. Webster, Kriging: a Method of Interpolation for Geographical Information System. Int. J. Geographical Information Systems. Vol. 4, No. 3 (1990), pp. 313-332.

34. A.A. Самарский, A.B. Гулин, Численные методы. M.: Наука, 1989.

35. В.Е. Bayer, An Optimum Method For Two-level Rendition of Continuous-tone Pictures. Proceedings IEEE, International Conference on Communications, 1973. Vol. 26, pp. 11-15.1. Предметный указатель

36. Глубина ячейки, 38 Граничные условияна боковой границе области, 109 на дне водоёма, 109 на поверхности водоёма, 108 непроницаемости, 33 первого рода, 35 периодические, 331. Дерево сетки, 27множество элементов, 22

37. Кластеризация ячеек, 16 Конвективная диффузия, 94 Концентрация примеси уравнение, 94

38. Метод конечных объёмов, 87 Модель

39. Бакли-Леверетта, 93 переноса многокомпонентной примеси, 107 Модификация сеток, 54

40. Объединение ячеек, 21, 28, 54 Ограничениена количество ячеек, 56 на размер ячеек, 47, 56 Окрестность ячейки, 21, 47 Операция на сетке, 54разрешённая, 56 Осаждаемое вещество, 101

41. Пассивная примесь, 92 Поддерево, 25элемента, 26 Потомок, 25, 57 косвенный, 24, 60 прямой, 25 Примеси, 99гранулометрический состав, 99 Пространство адаптивных сеток, 27, 42, 55 метрика, 54 производное, 68

42. Разбиение ячейки, 20, 54 Размеры ячейкиограничение на размеры, 21 Расстояниеструктурное, 541. Сеткадочерняя, 19 Сеточная функция, 21, 27 интеграл, 51фиксация функции, 671. Соседние ячейкиалгоритм нахождения, 39 малые, 30, 42 обычные, 29, 48

43. Сохранение данных сетки, 78

44. Структура четверичного дерева, 38 функция структуры, 22, 26

45. Сшивка адаптивных сеток, 32схемная вязкость, 73

46. Управляющий функционал алгоритма, 66

47. Уровень измельчения сетки, 151. Фотическая зона, 1051. Фракции, 99

48. Четверичное дерево, 19, 221. Элементкорневой, 23 листовой, 23 родительский, 231. Ячейкадочерняя, 20, 28 корневая, 20 листовая, 26 родительская, 20