Содержание к диссертации
Введение
Глава 1. Фильтрация и перенос в насыщенной пористой среде: основные понятия и математические модели 14
1.1. Модель фильтрации 14
1.1.1. Параметры пористых сред 14
1.1.2. Закон Дарси 15
1.1.3. Закон сохранения массы 16
1.1.4. Граничные и начальные условия 17
1.2. Модель переноса примесей в пористой среде 18
1.2.1. Конвекция 18
1.2.2. Диффузия и дисперсия 18
1.2.3. Уравнение конвекции-диффузии 19
1.2.4. Граничные и начальные условия 21
Глава 2. Численные методы решения задач переноса примесей в пористых средах 23
2.1. Тетраэдральная сетка и сеточные пространства 24
2.2. Семейство монотонных методов конечных объемов для численного решения диффузионных задач
2.2.1. Формулировка методов 25
2.2.2. Свойства методов 32
2.3. Методы дискретизации задач конвекции-диффузии 36
2.3.1. Схема Жаффре: РКЭ для оператора конвекции и МСКЭ для оператора диффузии 37
2.3.2. Схема расщепления по физическим процессам
2.4. Результаты численных экспериментов: задачи с гладким решением 45
2.4.1. Используемые обозначения 45
2.4.2. Задача с доминирующей диффузией 47
2.4.3. Задача с доминирующей конвекцией 49
2.4.4. Задача с полным тензором диффузии 51
2.5. Результаты численных экспериментов: задача о переносе фрон
та концентрации 52
Глава 3. Методы повышения эффективности вычислений 59
3.1. Параллелизация алгоритма 59
3.1.1. Разбиение области по процессорам 59
3.1.2. Параллелизация локальных шагов 59
3.1.3. Параллелизация итерационного метода решения глобальной системы 60
3.1.4. Результаты численного эксперимента 62
3.2. Применение динамических сеток 68
3.2.1. Технология перестроения сетки 68
3.2.2. Переинтерполяция решения с одной сетки на другую 70
3.2.3. Результаты численных экспериментов 71
3.3. Двухуровневый переобуславливатель для задач диффузии с анизотропным неоднородным тензором диффузии 75
3.3.1. Модельная задача, дискретизация и необходимые обозначения 75
3.3.2. Построение переобуславливателя 78
3.3.3. Реализация алгоритма 79
3.3.4. Результаты численных экспериментов
Глава 4. Решение прикладной трехмерной задачи о распростра нении ядерных загрязнений 86
4.1. Постановка задачи 86
4.2. Выбор расчетной тетраэдральной сетки 90
4.3. Выбор расчетных схем 93
4.4. Результаты расчета гидравлического напора 96
4.5. Результаты расчета распространения загрязнения 99
Заключение 107
Литература 108
- Граничные и начальные условия
- Схема Жаффре: РКЭ для оператора конвекции и МСКЭ для оператора диффузии
- Параллелизация итерационного метода решения глобальной системы
- Выбор расчетной тетраэдральной сетки
Введение к работе
Актуальность темы. Необходимость в решении задачи переноса примесей в пористых средах возникает при математическом моделировании распространения ядерных или химических загрязнений в геологических пластах Источниками таких загрязнений служат подземные хранилища радиоактивных отходов, сооружаемые в ряде стран с развитой ядерной энергетикой Интерес к этим задачам и их актуальность обусловлены необходимостью прогноза миграции загрязнений с целью обеспечения безопасности жизнедеятельности человека
Цель диссертационной работы. Целью данной диссертационной работы является разработка вычислительных схем и создание программного комплекса для трехмерного численного моделирования задач переноса примесей в пористых средах с реальными параметрами
Научная новизна. В работе предложено новое семейство методов конечных объемов для численного решения диффузионных задач с полным неоднородным тензором диффузии на неструктурированных тетраэдральных сетках, доказана теорема об их монотонности Создана технология решения нестационарных задач конвекции-диффузии на динамически перестраиваемых сетках Проведены следующие исследования сравнение схем расщепления для нестационарного уравнения конвекции-диффузии с традиционными схемами, тестирование нового переобуслав-ливателя для диффузионных задач и его сравнение с существующими Решена прикладная трехмерная задача с реальными коэффициентами
Достоверность результатов. Достоверность результатов диссертационной работы обоснована использованием теории численных методов и строгих математических выводов, результатами вычислительных экспериментов на тестовых примерах
Используемые методы. В исследовании применяются следующие методы Построение математических моделей рассматриваемых задач основывается на законах Дарси и сохранения массы Разработка нового семейства трехмерных методов конечных объемов строится на идее двумерного монотонного метода конечных объемов К Ле Потье Для дискретизации задач используются метод конечных элементов, метод конечных объемов, метод смешанных конечных элементов (в гибридной постановке), схемы расщепления оператора
Решение возникающих линейных систем осуществляется итерационными методами на подпространствах Крылова Для переобуславливания
систем применяются метод неполной факторизации ILU2 И Е Капорина и параллельный метод декомпозиции области Ю В Василевского Вывод и анализ нового переобуславливателя для диффузионных задач основан на теории двухуровневых переобуславливателей Ю А Кузнецова
При распараллеливании алгоритмов используется метод инерци-альной бисекции для разделения сетки между процессорами и MPI-библиотека для организации межпроцессорных обменов
Теоретическая и практическая ценность. Теоретическая ценность работы состоит в разработке нового семейства методов конечных объемов для дискретизации диффузионного оператора Метод обладает свойством монотонности в смысле неотрицательности получаемого численного решения стационарного или нестационарного уравнения диффузии с соответствующими начально-краевыми условиями и источником
Практическая ценность заключается в создании программного комплекса на языке Fortran для численного трехмерного моделирования распространения ядерных или химических загрязнений в геологических пластах Программа имеет параллельную реализацию, что важно при решении больших задач Решена практическая задача о распространении ядерных загрязнений из подземного хранилища радиоактивных отходов
Апробация работы Результаты диссертационной работы докладывались автором и обсуждались на научных семинарах Института вычислительной математики РАН, Института математического моделирования РАН, в университете г Линчопинг (Швеция) и на следующих конференциях
П-ая международная конференция по матричным методам и операторным уравнениям, Москва, июль 2007г
Конференция SIAM по математическим и вычислительным аспектам геологии, США, Санта-Фе, март 2007г
XLIX-ая научная конференция МФТИ "Современные проблемы фундаментальных и прикладных наук", Москва, декабрь 2006г
VII -ой Международный научный симпозиум по итерационным методам IMACS, США, Колледж Стейшн, университет Texas А&М, ноябрь 2006г
Всероссийская молодежная школа-конференция "Численные методы решения задач математической физики", Казань, июнь-июль 2006г
Международная конференция "Тихонов и современная математика", Москва, июнь 2006г
Конференция "Ломоносовские чтения", Москва, апрель 200бг
Международный симпозиум по проблемам моделирования захоронений ядерных отходов MOMAS, Франция, Марсель, ноябрь 2005г
Конференция "Тихоновские чтения" , Москва, октябрь 2005г
Публикации По теме диссертации опубликовано 5 работ 2 — в рецензируемых журналах (входят в перечень ВАК), 1 — в сборнике научных трудов, 2 — в материалах конференций Список работ приведен в конце автореферата
Личный вклад автора В совместной работе [2] вклад автора заключался в программной реализации переобуславливателя и постановке численных экспериментов, в работе [5] — в разработке параллельного алгоритма, его программной реализации и проведении расчетов на многопроцессорной ЭВМ
Структура работы. Диссертационная работа состоит из введения, четырех глав, заключения и списка литературы из 62 наименований, включает 38 рисунков и 21 таблицу Общий объем диссертации - 115 страниц
Граничные и начальные условия
Жидкость, занимающая поровое пространство в насыщенной пористой среде, способна переносить примеси в растворенном виде. Процесс конвекции заключается в распространении раствора примеси вместе с общим фильтрационным, или конвективным, потоком. Уравнение конвекции для насыщенной пористой среды выглядит следующим образом: ЯГ и, ± =-V (ОТ). (1.12) Здесь С — концентрация примеси в растворе,[ML-3]. Пористость в левой части уравнения появляется вследствие того, что конвективный перенос примеси происходит со скоростью потока жидкости в порах, которая получается делением скорости фильтрационного потока и на величину пористости и [5].
Молекулярная диффузия и кинематическая дисперсия отражают процессы смешивания жидкостей, заполняющих поры. Эффект молекулярной диффузии связан с наличием грандиента концентрации примеси в жидкости. Этот процесс не зависит от скорости фильтрационного потока и приводит к перераспределению примеси из зон высокой концентрации в зоны низкой концентрации.
Кинематическая дисперсия заключается в перемешивании жидкостей вследствие неоднородности поля фильтрационного потока на микроскопическом уровне. Уравнение диффузии-дисперсии в пористой среде выглядит следующим образом: дС u— = V-VVC, (1.13) где V, [L2T l], — симметричная положительно определенная матрица размера 3 х 3, в дальнейшем называемая тензором диффузии и состоящая из двух компонентов: V = VC + VM. (1.14) Тензоры Т с и Т м отвечают за кинематическую дисперсию и молекулярную диффузию, соответственно. Т м — диагональная матрица 3x3. Для тензора Т с обычно используется формула [38],[29]: Vc = \\й\\(аіЕ(и) + at(I - 5(iZ))). (1.15) В (1.15) введены следующие обозначения: I — единичная матрица размера 3x3; Е(и) — тензор с компонентами Eij(u) — w, hj Є {1,2,3}; а/ — коэффициент продольной дисперсии; at — коэффициент поперечной дисперсии.
В отсутствие химических реакций, радиоактивного распада и взаимодействия примеси с твердым скелетом масса примеси остается постоянной и ее перенос описывается уравнением конвекции-диффузии в пористой среде: ВС и— = - V (Си) + V VVC. (1.16) at Перенос примеси может замедляться под действием эффекта адсорбции [5]: молекулы растворенной в жидкости примеси временно закрепляются (адсорбируются) на поверхности пор. Между раствором и поверхностью устанавливается динамическое равновесие, и доля адсорбированной примеси в твердом скелете са (безразмерная величина) зависит от его концентрации в растворе. Зависимость са(С) носит название изотермы адсорбции и часто предполагается линейной: са = коС, ко — коэффициент распределения, [M_1L3]. Общая масса примеси в единичном объеме определяется выражением: шС + р.{\ -ш)са =w(l + р3]- -кЛ С = uRC, (1.17) где ps — плотность твердого вещества в пористой среде; R = (l + Ps f-ho) — коэффициент задержки. Уравнение переноса примеси в насыщенной пористой среде с учетом адсорбции записывается следующим образом: дС uR— = -V-{Cu) + 4-VVC. (1.18) Очевидно, что при R 1 процессы конвекции-диффузии замедляются эффектом адсорбции.
Для расчета задач распространения ядерных загрязнений требуется учет радиоактивного распада веществ. Пусть Тд — период полураспада радиоактивного загрязнителя, тогда при наличии источника загрязнения /, [МЬ гТ 1] закон сохранения массы примеси выглядит следующим образом: uR ( + \с\ = -V (Си) + V VVC + /, (1.19) где Л = In 2/Td — постоянная радиоактивного распада. 1.2.4. Граничные и начальные условия
Для нестационарного уравнения сохранения массы вида (1.19) требуется задание начальной концентрации при t — 0 и граничных условий. В задачах переноса примесей [29] обычно используются следующие типы граничных условий: Граничные условия типа Дирихле, определяющие концентрацию на гра нице То C(x,t) = CD(x,t), xeTD. Граничные условия типа Неймана, определяющие диффузионный по ток через границу Гдг . -V = rN(x,t),xeTN. Использование такого типа граничных условий возможно на тех участках границы, где и п 0, п — внешняя нормаль к границе. Обычно эти условия ставятся на удаленной границе, где считается, что градиент решения близок к нулю, —Т Щ — 0. Общий поток через границу Г А/: дС -V— + Сй-п = тм(х, t), х Є Гм. on Такие условия с нулевой правой частью характерны для непроницаемых границ. Выводы по главе 1 При моделировании переноса примесей в насыщенных пористых средах необходимо решение двух задач. Во-первых, на каждом шаге дискретизации по времени нужно определить фильтрационные потоков с помощью закона Дарси (1.3) и уравнения неразрывности (1.6). Во-вторых, требуется решать задачу конвекции-диффузии типа (1.16), возможно, с дополнительными членами радиоактивного распада или химической реакции (1.19). В пренебрежении изменением плотности жидкой фазы, заполняющей поры, эти задачи решаются независимо друг от друга, сначала находятся фильтрационные потоки, затем они используются для расчета переноса примеси.
Схема Жаффре: РКЭ для оператора конвекции и МСКЭ для оператора диффузии
Исключение неизвестных производится поэлементно, на каждом тетраэдре. В результате получена система с положительно определенной симметричной матрицей порядка Njrr Для ее решения используется метод сопряженных градиентов с переобуславливателем ILU2 [43] (модифицированный метод неполной факторизации второго порядка с двумя порогами). Определив из (2.35) А"4"1, восстановим средние концентрации и потоки локально по формулам с = -(BrA-1B + G)-1(F2-BTA-1F1+BrA-1CA"+1), (2.36а) r„+1 = A-i(Fl _ Вс _ сдп+ij (2.36Ь)
Важной особенностью представленной схемы является то, что расщепление на физические процессы в ней появляется только в процессе разложения пространства тестовых функций на ортогональные подпространства при решении системы уравнений (2.28Ь)-(2.28с). На каждом подшаге схемы (2.28) уравнение конвекции-диффузии аппроксимируется полностью. В силу этого, в метод не вносится дополнительная аппроксимационная ошибка, что позволяет рассматривать его формально как метод второго порядка по времени и по пространству. Противопотоковая аппроксимация конвективного члена также не снижает порядка аппроксимации благодаря использованию разрывных конечных элементов.
Заметим, что вышеописанное разложение пространства Wh(T) невозможно на тетраэдрах, имеющих точки на границе с условиями типа Дирихле. В этом случае шаги (2.31а)-(2.31с) выполняются формально, при этом граница области разделяется на две части: сЮ = dQ,mUdQ,ut, где дО,гп = {є;, і Є Тв и пЄі 0} и dflout = {ЄІ, І Є ТВ и- пЄі 0}, причем пЄі является внешней нормалью к dQ. На части границы 80, граничные условия выполняются строго, а на dQ0Ut граничные условия выполняются только в смысле средних концентраций на гранях Л/г Є A/t. Численные эксперименты подтверждают эффективность и устойчивость такой трактовки граничных условий.
Как и в предыдущей схеме, концентрация приближается кусочно-линейными разрывными функциями из пространства Wh. Схема использует расщепление по физическим процессам, операторы диффузии и конвекции разнесены по подшагам, и на каждом подшаге решается неполное уравнение [23]:
Для аппроксимации конвективного оператора выбрана явная схема предиктор-корректор с регуляризацией против потока в корректоре. В предикторе (2.37а) вычисляется промежуточная концентрация Ch 2, в корректоре — концентрация С а на основе конвективных потоков на промежуточном слое. В интеграле по границе след Ch Лц4 берется на тетраэдре, лежащем выше по потоку. К С а применяется процедура ограничения наклона (2.37с), аналогичная (2.28d), и далее с помощью неявной схемы (2.37(1) вычисляется добавка к средней концентрации за счет диффузионных потоков г 1 через грани тетраэдра Т(г — 1,...,4). Для нахождения C+i и г 1 на диффузионном шаге (2.37d) может быть применен новый нелинейный МКО, сформулированный в 2.2, или МСКЭ. Поясним подробнее использование МКО с коэффициентами (2.12) (т.е., параметр ре принят равным нулю для всех граней є Є dh).
Для реализации шага (2.37d) найдем проекцию с решения С+ а на точки множества В (определено в 2.2) и решим нелинейную задачу МКО для с - , искомых концентраций в точках В: (V + A(cdiff)At) cdiff = Vc. (2.38) Здесь V — диагональная матрица объемов элементов, А(с?1Н) - несимметричная матрица МКО, элементы которой зависят от с Л Следствие 1 2.2 гарантирует неотрицательность решения сА , если с неотрицательно. После вычисления cdl воспользуемся формулой и найдем добавку к средним концентрациям за счет диффузионных потоков, требуемую в (2.37е).
Важным достоинством данной схемы по отношению к схеме Жаффре (см. раздел 2.3.1) является возможность использования разных шагов по времени для процессов конвекции и диффузии, аналогично [33], [23]. Более подробно это описано в 4.3. Шаг по времени для конвективных процессов определя ется требованиями устойчивости схемы, в то время как для диффузионных процессов — только требованиями точности. При равных шагах для конвекции и диффузии (как и в схеме Ж.Жаффре) реализация (2.37d) занимает около 80% времени, и увеличение шага для диффузии дает существенную экономию времени счета. Достоинством схемы Жаффре остается формально второй порядок аппроксимации, так как расщепление в ней появляется вследствие разложения пространства тестовых функций, и на каждом шаге решается полное уравнение.
Параллелизация итерационного метода решения глобальной системы
Персобуславливатель Li определяется неявно через уравнение L2 = TLTT, (3.5) где L — последовательный ILU2 или ILUt переобуславливатель для М. Отметим, что Li — сглаживатель в пространстве агрегированных векторов V — {г R": v = Tv, v Є R"}, которые могут целиком обрабатываться на одном процессоре, поскольку п п. В силу того, что персобуславливатель L несимметричен, будем использовать метод минимальных невязок для решения переобусловленной системы (12). Анализ переобуславливателя L представлен в работе [60], детальное рассмотрение одной его параллельной реализации можно найти в [54]. Параллельные обмены при применении рассматриваемого переобуславливателя сводятся к: локальной пересылке между процессорами-соседями для распределения результата применения L\\ локальной пересылке при вычислении невязки, глобальной пересылке на корневой процессор агрегированного вектора для применения L.
Постановка задачи Couplexl. Задание Couplexl [29] было предложено французским агентством по радиоактивным отходам Андра (www.andra.fr). Изначально задача формулировалась как двумерная, мы же расширяем ее до трех измерений. Геометрия расчетной области представлена на рис.3.1. Она представляет собой параллелепипед Q = (0,25000) х (0, 25) х (0,695) в метрах. Геологические слои расположены следующим образом:
Источник загрязнений (контейнер с радиактивными отходами) считается прямоугольным параллелепипедом, расположенным в глине R = {(х, у, z) е (18440,21680) х (0,25) х (244,250)}. Фильтрационные потоки подчиняются законам Дарси (1.3) и неразрывности в виде (1.8). Они посчитаны в рамках отдельной задачи с помощью МСКЭ (см. 4.3). Распределение гидравлического напора показано на рис. 3.2. Оно не зависит от координаты у, линии уровня идут по возрастанию слева направо. Скорость фильтрации достаточно велика только в известняке и юрском известняке. В других слоях распространение загрязнения определяется преимущественно диффузионным процессом. Подробное описание поля скоростей для задачи Couplexl можно найти в [39].
В качестве радиоактивного загрязнителя рассматривается Йод 129 (1291), мощность источника задана во времени (рис. 3.3) и распределяется равномерно по площади вертикального сечения источника (3240 х 6м2), постоян "і і 1 1 1 1 і 1 1 1 1 г на вдоль оси Y. Распространение загрязнения описывается уравнением конвекции-диффузии-распада (1.19). Период полураспада 1291 составляет 1.57 107лет. Тензор диффузии-дисперсии V определяется из соотношений (1-14), (1.15).
Значения коэффициентов, входящих в (1.19), представлены в таблице 3.1. За начальное время t = 0 принят момент, когда контейнер начинает течь, 100000 Мощность источника как функция времени, (кг/год). что задает нулевые начальные условия для концентрации во всей области Q. Граничные условия не зависят от времени: дС — = 0 на {0} х (0,25) х (295, 595), {0}х (0,25) х (0,200), (0,25000) х {0,25} х (0,695). VVC-п-Сй-п = 0 на (0,25000) х (0,25) х {0}. С — 0 на остальных границах. В задаче требуется проиллюстрировать поверхности уровня концентрации Ю-12, Ю-10,10"8,10 6кг/м3 в разные моменты времени.
Результаты расчетов. Для расчета задачи Couplexl была создана тет-раэдальная сетка, содержащая 786 тысяч тетраэдров и 1.65 миллиона граней. Сетка имеет сгущения к источнику загрязнений и границе между глиной и известняком. Шаг схемы по времени составляет 25 лет (первые 10 шагов — по 100 лет), что обеспечивает устойчивость схемы (2.28) па данной пространственной сетке. В каждом эксперименте делается 1000 шагов. Сравниваются времена счета на разном числе процессоров, с разными вариантами переобу-славливания в подобластях.
Расчеты проводились на суперкомпьютере MVS-15000, принадлежащем Межведомственному Суперкомпыотерному Центру РАН. MVS-15000 имеет 864 вычислительных узла. Каждый узел содержит два 64-разрядных процессора IBM PowerPC 970 с тактовой частотой 1600МГц и 4Гб оперативной памяти. Для обменов данными при счете задач используется сеть Myrinet с производительностью около 450Мбайт/с в режиме обмена "точка-точка".
Результаты расчетов приведены в таблице 3.2. На рисунках 3.4, 3.5 представлено решение задачи в моменты времени 13250 лет (500 шагов) и 25750 лет (1000 шагов). На них изображены поверхности уровня, соответствующие концентрации Ю-12, Ю-10, Ю-8,10_6кг/м3, черным цветом закрашен источник загрязнения. Из рисунков видно, что в глине загрязнение распространяется в вертикальном направлении. Доходя до известняка, оно сносится сильным горизонтальным конвективным потоком. В ячейках таблицы 3.2 представлены время расчета на разном числе процессоров (Т4,Г8,Тіб), среднее число итераций на шаг по времени (Niter) и показатели эффективности па-раллелизации Si (і Є {8,16}), определяемые по формуле Si = A » = 8,16 . (3.6) ±i %
Критерием остановки итерационного метода является относительное падение евклидовой нормы невязки на 6 порядков. Малый рост числа итераций при увеличении числа процессоров говорит об эффективности параллельного пе-реобуславливателя. Время расчета при увеличении числа процессоров уменьшается, однако эффективность параллелизации падает, так как повышается доля обменов данными в общем объеме работы. В то же время отметим, что с увеличением числа процессоров растет объем доступной памяти, что позволяет выбрать лучшие параметры переобуславливателя и значительно ускорить
Выбор расчетной тетраэдральной сетки
Решение при Т = 50000 лет (слева) и Т = 100000 лет (справа) в сечении плоскостью Z = 35м. показан черными стрелками на рис.4.14 (справа). Динамика других изопо-верхностей концентрации аналогична, к примеру, уровень С = 10 12кг/м3 выходит в Оксфорд Нр1-Нр4 приблизительно в 15000 лет и достигает непокрытого киммериджа в 60000 лет; для уровня С = 10 10кг/м3 эти времена составляют 45000 лет и 180000 лет соответственно.
Влияние фильтрационных потоков в пластах Оксфорд СЗа-СЗЬ, L2a-L2b, Hpl-Hp4 показано на рис.4.15, 4.16. На них представлены изолинии концентрации загрязнения, соответствующие 10 14, Ю-12,10 10кг/м3 в сечении области плоскостью Z = 35м.
Доходящие до верхних непокрытых слоев загрязнения оказываются на поверхности области в течение нескольких лет. Эти процессы проиллюстрированы на рис.4.17, где показано поле концентрации в непокрытом киммеридже и титоне, получаемое в результате моделирования с нулевыми начальными условиями и стационарными граничными условиями первого рода на T F, соответствующими решению в нижних слоях на этой границе при Т = 350000 лет. Рис.4.17 соответствует решению в непокрытых пластах через 5 лет с
Распределение времени по процедурам при параллельном расчете модели медленного распространения.
На всех рисунках решение представлено в логарифмическом масштабе. Полученные результаты носят характер прогноза, поскольку численная диффузия, присутствующая в схеме, способна привести к завышению миграционной способности радионуклидов и ускорению процессов переноса по отношению к реальным.
Решение задачи в рамках модели медленного переноса заняло 29 часов при использовании 12 процессоров на многопроцессорной ЭВМ MVS-50000 Межведомственного суперкомпьютерного центра. Каждый узел этого компьютера имеет 4 процессора Intel Itanium II и 4Гб разделяемой памяти. В таблице 4.3 показано распределение расчетного времени по подшагам схемы (4.4). При решении линейной системы на шаге (4.4d) критерием остановки итерационного метода было относительное падение евклидовой нормы невязки на 6 порядков. Для этого требовалось в среднем 14 итераций метода бисопряжен-ных градиентов с переобуславливателем. Большое время, затрачиваемое на расчет диффузионного шага (4.4d), оправдывает использование схемы расщепления с разными шагами по времени для диффузии и конвекции. Такое расщепление позволяет на порядок сократить время расчета относительно схемы с одинаковыми шагами. Для модели быстрого переноса использовалась сеть уровня Gigabit Ethernet из трех двухъядерных компьютеров Intel Pentium D (ИВМ РАН), время расчета составило порядка 10 минут.
В главе проведено численное моделирование прикладной задачи Andra 3D, включающее построение тетраэдральной сетки и раздельное решение задач фильтрации и конвекции-диффузии. Рассмотрены три подхода к построению тетраэдральных сеток, выбран наиболее адекватный метод для данной задачи и создана сетка, содержащая 2.5 миллиона ячеек. Найдено поле фильтрационных потоков и гидравлический напор в условиях сильной анизотропии и неоднородности тензора фильтрации и расчетной сетки.
Решение задачи переноса также осложнялось сжатием сетки по вертикали, анизотропией и неоднородностью тензора диффузии, разными временными масштабами процессов в верхних и нижних геологических пластах. В результате расчетов получен прогноз времени и места выхода на земную поверхность загрязненных вод с концентрациями радиоактивного Йода 129, равными Ю-14, Ю-12, 10_10кг/м3 . Показана динамика распространения загрязнения в трехмерной области с течением времени, вплоть до 350000 лет.
Диссертационная работа посвящена разработке и исследованию численных методов, применяемых для моделирования трехмерных задач переноса примесей в пористых средах с полными и неоднородными тензорами диффузии и проницаемости. Для дискретизации расчетных областей используются тетраэдральные неструктурированные сетки. Главным результатом работы стало создание программного комплекса для численного решения задач переноса примесей в пористых средах и его применение для моделирования прикладной задачи.
Основная теоретическая ценность работы заключается в разработке нового метода конечных объемов для дискретизации диффузионного оператора с полным анизотропным неоднородным тензором диффузии. Метод обладает свойством монотонности в смысле диагонального преобладания по столбцам в матрице аппроксимации, что гарантирует неотрицательность решения при неотрицательных краевых условиях первого рода, неположительных — второго рода и неотрицательных источниках. При однородных краевых условиях типа Дирихле на всей границе области и неотрицательных источниках метод отвечает дискретному принципу максимума.
Практическая ценность работы состоит в создании программного комплекса на языке Fortran для численного трехмерного моделирования распространения химических или ядерных загрязнений в геологических пластах. В него включены блоки построения тетраэдральной сетки, решения фильтрационной задачи и моделирования нестационарной задачи переноса примеси в пористой среде. Реализована параллельная версия двух последних блоков, что важно при решении больших задач с миллионами неизвестных. Решена прикладная задача о распространении ядерных загрязнений из хранилища радиоактивных отходов.
