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

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

Капырин Иван Викторович

Трехмерное моделирование процессов переноса примесей в пористых средах сложной структуры

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

АВТОРЕФЕРАТ диссертации на соискание ученой степени кандидата физико-математических наук

Москва 2007

003161717

Работа выполнена в Институте вычислительной математики РАН

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

Ю В Василевский

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

В И Агошков

Кандидат технических наук, А В Расторгуев

Ведущая организация. Институт математического

моделирования РАН

Защита состоится 9 ноября 2007 года в 15 часов на заседании диссертационного совета Д 002 045 01 в Институте вычислительной математики РАН по адресу 119333, Москва, ул Губкина, д 8, ауд 727

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

Автореферат разослан октября 2007 года

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

диссертационного совета

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

ГА Бочаров

ОБЩАЯ ХАРАКТЕРИСТИКА РАБОТЫ

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

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

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

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

Используемые методы. В исследовании применяются следующие методы Построение математических моделей рассматриваемых задач основывается на законах Дарси и сохранения массы Разработка нового семейства трехмерных методов конечных объемов строится на идее двумерного монотонного метода конечных объемов К Ле Потье Для дискретизации задач используются метод конечных элементов, метод конечных объемов, метод смешанных конечных элементов (в гибридной постановке), схемы расщепления оператора

Решение возникающих линейных систем осуществляется итерационными методами на подпространствах Крылова Для переобуславливания

систем применяются метод неполной факторизации 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 страниц

СОДЕРЖАНИЕ РАБОТЫ

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

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

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

Модель фильтрации основана на законе Дарси и законе сохранения массы жидкости Закон Дарси устанавливает зависимость фильтрационных потоков от градиента гидравлического напора

й=-^Н = --{\7р + рдЧг), (1)

I1

где и — вектор скорости фильтрационного потока, р — плотность жидкости, которая предполагается постоянной в среде, д — гравитационная постоянная, р — давление, К = — коэффициент фильтрации, зависящий от свойств пористой среды и жидкости, г — геометрический напор (высота точки), Н = + г — гидравлический напор, р, — динамическая вязкость жидкости Закон сохранения массы жидкости постулирует, что изменение массы жидкости в элементарном объеме определяется изменением за счет сжимаемости и за счет внешних источников или стоков с массовой интенсивностью /

-V = (2)

Система (1),(2) дополняется граничными условиями, которые могут быть типа Дирихле, типа Неймана и смешанными

Во втором, параграфе описана модель переноса примеси в пористой среде Последовательно дано описание физических явлений конвекции,

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

где Д — коэффициент задержки, и> — пористость, С — концентрация примеси в растворе, Л — постоянная радиораспада, Т> — тензор диффузии-дисперсии (симметричная положительно определенная матрица 3x3) На границах области возможна постановка граничных условий типа Дирихле, типа Неймана или задание общего потока примеси через границу В пренебрежении изменением плотности жидкой фазы, заполняющей поры, задачи фильтрации и переноса решаются независимо друг от друга сначала находятся фильтрационные потоки, затем они используются для расчета распространения примеси

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

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

Во втором параграфе дана формулировка нового метода конечных

(3)

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

V г = /, (4)

г = -DVC в Q, (5)

С\То = д°(х), ? Й|Глг = gN{x)

Здесь г— диффузионный поток, fi G R3 — расчетная область с границей дй = Гд U Глг

Построим в области fi сетку е^, содержащую Nt тетраэдров Для каждого тетраэдра Т в сетке £h вводится точка Хт — носитель степени свободы на Т Пусть Т = ABCD с гранями a, b, с, d, противолежащими A,B,C,D, соответственно, Ra,Rb,Rc>Rd ~ радиус-векторы соответствующих вершин Т, а векторы ña,ñb,ñc,ñd — внешние нормали к граням, длины нормалей численно равны площадям сответствующих граней Для радиус-вектора Вх1 точки Хт положим

5 = ДаЦтЦр + Йв\\щ\\у + Rc\\nc\\v + Дд1Ы|р Xl \\na\\v+\\ñb\\v+\\ñc\\v + \\ñd\\v ' W

где \\ñp\\v = y/{T>ñ0,ñ0), /3 e {a, b, c, d}

Интегрируя по каждому T E £h закон сохранения массы (4) с использованием формулы Гаусса-Остроградского и переходя к дискретной формулировке, получаем

Е

еедТ т

Ге пе = J f áx VT € eh, (7)

где пе — внешняя нормаль к грани е тетраэдра Т, длина которой численно равна площади соответствующей грани, те \пе\ — |е|, а ге пе ~ аппроксимация диффузионного потока через грань е Главной идеей метода является способ построения аппроксимации диффузионного потока ге пе Для этого рассматриваются два соседних тетраэдра Т+ = АО^О^О?, и Т_ = ВО1О2О3 (см рис 1) Точки Х+,Х_ — соответствующие Т+ и Т_ точки-носители, М —центр масс общей грани е, е = О1О2О3 С помощью интегрирования (5) по тетраэдрам вида Х+МО0О^ Х^МО^Ок строятся аппроксимации диффузионных потоков гг пе/3, % € {1,2,3}, через каждый из треугольников МО^О^ {г,3,к — различны) Общий поток через грань представляется линейной комбинацией этих трех потоков

Ге Пе = (4?! пе + Ц2 Г'2 пе + Пе (8)

о,

Рис 1 Геометрические построения для нелинейного метода конечных объемов

Он зависит от значений концентрации Сх+, Сх_ в точках-носителях Х+,Х- и Со, в точках Ог, г € {1,2,3} Значения Со, с помощью линейной интерполяции выражаются через концентрацию в точках-носителях

Для определения коэффициентов $, г = 1,2,3, накладываются ограничения на диффузионный поток через грань е (8)

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

з

• Шаблон аппроксимации потока является двухточечным

Требования двухточечной аппроксимации потока определяют семейство решений (ре), г е {1,2,3}, с параметром ре

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

(9)

А(СХ)СХ =

(10)

где Сх —вектор неизвестных концентраций в точках-носителях Хт, размера Nt Матрица А(Сх) — разреженная, имеющая не более 5 ненулевых элементов в каждой строке

Система (10) решается с помощью итерационного алгоритма Пикара

Л(Скх)Ск/1 = F (11)

с некоторым начальным приближением

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

Теорема. Положим в задаче (4) правую часть f(x) > 0 в ft, граничные условия gD{x) > 0 на Го, gN(x) < 0 на Гn Пусть (10) — соответствующая задаче (4) нелинейная система МКО, для которой точки-носитель степеней свободы на тетраэдрах определяются формулой (6), начальное приближение > 0, и на, каждой итерации метода Пи-

кара (11) для любой внутренней грани е выбираются неотрицательные значения ¿4®, % 6 {1,2,3} (последнее всегда возможно) Тогда все итерационные приближения к Сх неотрицательны,

{Ckx\> 0, i = l, , Nt, Vfc> 0

Доказательство теоремы основано на том факте, что в условиях теоремы матрица А{Сх) монотонна, а правая часть F неотрицательна

В конце параграфа приведены примеры построения на основе нового МКО схем для нестационарного уравнения диффузии Схемы сохраняют неотрицательность решения

В третьем параграфе главы рассматриваются две схемы расщепления для аппроксимации нестационарного уравнения конвекции-диффузии схема Ж Жаффре, использующая разрывные конечные элементы (РКЭ) для оператора конвекции и метод смешанных конечных элементов (МСКЭ) для оператора диффузии, и новая схема расщепления с монотонной МКО-аппроксимацией оператора диффузии и РКЭ для оператора конвекции

В схеме Жаффре аппроксимация по времени построена следующим образом используется явная схема предиктор-корректор для конвективного оператора и схема Кранка-Николсон — для диффузионного Для стабилизации схемы в корректоре применяется противопотоковая аппроксимация конвективного члена Использование РКЭ требует на каждом шаге по времени модифицировать решение с помощью ограничителя

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

Во второй схеме реализовано расщепление по физическим компонентам, операторы диффузии и конвекции разнесены по подшагам, и на каждом подшаге решается неполное уравнение Для аппроксимации конвективного члена снова применяется метод РКЭ и явная схема предиктор-корректор, для диффузионного члена — новый монотонный МКО (или МСКЭ) и неявная схема Важным достоинством данной схемы по отношению к схеме Жаффре является возможность использования разных шагов по времени для процессов конвекции и диффузии Шаг по времени для конвективных процессов определяется требованиями устойчивости схемы, в то время как для диффузионных процессов — только требованиями точности При равных шагах для конвекции и диффузии (как в схеме Ж Жаффре) реализация диффузионного под-шага занимает около 80% времени, и увеличение шага для диффузии дает существенную экономию времени счета

В четвертом и пятом параграфах проводится сравнение вышеописанных схем расщепления с двумя неявными аппроксимациями уравнения конвекции-диффузии неявной трехслойной схемой второго порядка по времени (ВБР) для МКЭ первого порядка с БИРв стабилизацией и ВБР-схемой МСКЭ с противопотоковой аппроксимацией конвективного члена Четвертый параграф посвящен тестированию на задачах с гладким решением, а пят,ый — на задаче о переносе фронта концентрации Схемы расщепления показывают порядок аппроксимации концентрации, близкий ко второму, и первый порядок аппроксимации потоков Использование разрывных конечных элементов делает их наиболее приспособленными к задачам переноса фронта концентрации, в которых они показывают низкую численную диффузию Схема Жаффре в применении к задаче с негладким решением демонстрирует очень малые нарушения монотонности Новая схема с нелинейным МКО дает результаты, близ-

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

Третья глава посвящена методам повышения эффективности вычислений В первом параграфе описан алгоритм распараллеливания расчетов применительно к схеме Жаффре Основные этапы параллели-зации следующие разбиение области по процессорам с помощью метода инерциальной бисекции, параллелизация локальных шагов, па-раллелизация итерационного метода решения СЛАУ Для переобуслав-ливания линейных систем используется метод декомпозиции области Ю В Василевского Проведены численные эксперименты на 4, 8, 16 процессорах, показывающие эффективность распараллеливания (время счета на 16 процессорах в 2 5-3 раза меньше, чем на 4 процессорах)

Во втором параграфе представлен алгоритм решения задач конвекции-диффузии на динамически перестраиваемых сетках Перестроение сетки осуществляется с помощью метода отмеченного ребра, позволяющего автоматически измельчать и разгрублять расчетную сетку в требуемых подобластях Предложены новые алгоритмы переинтерполяции решения с одной сетки на другую при измельчении или загрублении сетки Они сохраняют не только массу вещества, но и диффузионные потоки через грани сетки Метод протестирован на двух задачах переносе фронта от границы области и вращении непрерывного "всплеска" по кругу Показано, что применение динамических сеток позволяет существенно (в 4 раза для первой задачи, на 30% для второй) сократить время счета и понизить численную диффузию по сравнению с решением на статической сетке

В т,рет,ъем параграфе исследуется переобуславливатель для диффузионных задач с неоднородным анизотропным тензором диффузии, разработанный Ю А Кузнецовым в сотрудничестве с автором диссертации Переобуславливатель сформулирован для случая пластовых сред и матрицы дискретизации простейшего МКО на сетке Вороного Проведенные численные эксперименты подтверждают теоретический результат Ю А Кузнецова о независимости числа обусловленности переобусловленной матрицы МКО-аппроксимации от коэффициентов диффузии и шага сетки по вертикальной оси Сравнение переобуславливателя с двумя другими (алгебраическим многосеточным методом Штюбена и блочным методом Якоби) показывает его эффективность

В четвертой главе приведены результаты численного моделирования тестового задания Апска ЗБ, предложенного французским агент-

ством по ядерным отходам (www andra fr) Данные задачи соответствуют измерениям в зоне реального захоронения ядерных отходов

В первом параграфе главы дано описание геометрии и пластовой структуры расчетной области, определены физические параметры, необходимые для моделирования фильтрации и переноса, а также граничные условия задачи Задача обладает рядом трудностей анизотропия (2 порядка) и неоднородность (до 10 порядков) тензора фильтрации, сильная анизотропия области (диаметр области — 55км, глубина — 600м), неоднородность и анизотропия полного тензора диффузии, выклинивания геологических пластов, необходимость прогноза на большой период времени (сотни тысяч лет)

Второй параграф посвящен проблеме выбора расчетной сетки Рассмотрено три подхода к построению расчетных сеток Первый подход заключается в построении полностью неструктурированной, сгущающейся к источнику загрязнения сетки с помощью библиотеки ani3D, разработанной Ю В Василевским, К Н Липниковым и др Второй подход позволяет строить тетраэдральные сетки путем разбиения вертикальных треугольных призм Основания призм являются элементами двумерной треугольной сетки в проекции расчетной области на плоскость XY Третий подход, разработанный автором в сотрудничестве с Ю В Василевским и В Н Чугуновым, заключается в построении гексаэдральной сетки, разбиении ее элементов на тетраэдры и локальном измельчении в области источника загрязнений Недостатком первых двух подходов является наличие в сетке тупых двугранных углов, что приводит к крайне плохой обусловленности матриц аппроксимации К тому же, первый подход имеет существенные ограничения на число тетраэдров в сетке Для расчетов с помощью третьего подхода была построена сетка, содержащая 2 5 млн тетраэдров (около 5 млн неизвестных в возникающих линейных системах)

В третьем параграфе обоснован выбор схем для расчета задач фильтрации и переноса Для расчета фильтрационных потоков применяется схема гибридного МСКЭ При моделировании распространения загрязнения используется схема расщепления по физическим компонентам (из гл 2) с разными шагами по времени для конвективного и диффузионного процессов Для аппроксимации конвективного переноса в ней использованы РКЭ, для диффузионного — МСКЭ

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

ляется то, что скорость фильтрационных потоков в верхних, непокрытых, пластах значительно (на 4-5 порядков) выше, чем в нижних Это утверждение справедливо и для диффузионных потоков, так как тензор дисперсии линейно зависит от скорости фильтрации Расчетная область условно разделяется на подобласти состоящую из верхних пластов, и Ду, включающую все остальные пласты Разделение на подобласти используется при моделировании переноса загрязнений

В пятом параграфе приведены результаты моделирования распространения ядерных загрязнений (Йод 129) на период в 350000 лет Расчет переноса осложнялся сжатием сетки по вертикали, анизотропией и неоднородностью тензора диффузии, разными временными масштабами процессов в верхних и нижних геологических пластах Задача решается в рамках двух моделей модели медленного распространения в подобласти и модели быстрого распространения в подобласти Ор Показана динамика распространения изоповерхности концентрации С = 10~14кг/м3 вплоть до 350000 лет В результате расчетов получен прогноз времени и места выхода на земную поверхность загрязненных вод с концентрациями радиоактивного Йода 129, равными Ю-14, Ю-12, 10_10кг/м3 Вычисления проводились на параллельной многопроцессорной ЭВМ с использованием 12 процессоров в течение 29 часов Полученные данные о распределении времени счета по шагам расчетной схемы оправдывают использование разных шагов по времени для моделирования конвективного и диффузионного процессов

В заключении кратко сформулированы результаты диссертационной работы и сделаны выводы о ее теоретической и практической ценности

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

1 Разработано новое семейство монотонных методов конечных объемов для аппроксимации диффузионных задач с полным и неоднородным тензором диффузии на неструктурированных тетраэдральных сетках

2 Создан программный комплекс для решения трехмерных задач переноса примесей в пористых средах, включающий алгоритмы построения тетраэдральной сетки, расчета фильтрационных потоков и расчета переноса примеси (два последних алгоритма параллели-

зованы)

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

4 Проведено численное моделирование трехмерной задачи с реальными параметрами

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

1 Капырин И В Семейство монотонных методов численного решения трехмерных задач диффузии на неструктурированных тетраэдральных сетках // Доклады Академии Наук 2007 Т 614, №5, С 588-593

2 Kuznetsov Yu А , Boiarkme О V , Kapyrm I V , Yavich N В Numerical anahsys of a two-level preconditioner for the diffusion equation with an anisotropic diffusion tensor // Russian J Numer Anal Math Modelling 2007 Vol 22, №4, P 377-391

3 Капырин И В Об использовании динамических сеток при решении задач конвекции-диффузии // Труды математического центра им Н И Лобачевского 2006 Т 31, С 166-175

4 Капырин И В Численное решение задач переноса загрязнений в подземных водоносных слоях на тетраэдральных сетках с использованием многопроцессорных ЭВМ // Book of abstracts of international conference "Tikhonov and contemporary mathematics", section 2 Mathematical modelling M Изд-во МГУ 2006 С 91-92

5 Василевский Ю В , Капырин И В Параллельное трехмерное моделирование распространения примесей в пористых средах // Матричные методы и технологии решения больших задач М ИВМ РАН 2005 С 33-50

Заказ N° 42/10/07 Подписано в печать 02 10 2007 Тираж 75 экз Усд п л 0,75

ООО "Цифровичок", тел (495) 797-75-76, (495) 778-22-20 ■шгм с/г ги, е-тай т/о@с/г ги

Введение

Глава 1. Фильтрация и перенос в насыщенной пористой среде: основные понятия и математические модели.

1.1. Модель фильтрации.

1.1.1. Параметры пористых сред.

1.1.2. Закон Дарси.

1.1.3. Закон сохранения массы

1.1.4. Граничные и начальные условия.

1.2. Модель переноса примесей в пористой среде.

1.2.1. Конвекция.

1.2.2. Диффузия и дисперсия.

1.2.3. Уравнение конвекции-диффузии.

1.2.4. Граничные и начальные условия.

Глава 2. Численные методы решения задач переноса примесей в пористых средах

2.1. Тетраэдральная сетка и сеточные пространства

2.2. Семейство монотонных методов конечных объемов для численного решения диффузионных задач.

2.2.1. Формулировка методов.

2.2.2. Свойства методов.

2.3. Методы дискретизации задач конвекции-диффузии.

2.3.1. Схема Жаффре: РКЭ для оператора конвекции и МСКЭ для оператора диффузии.

2.3.2. Схема расщепления по физическим процессам.

2.4. Результаты численных экспериментов: задачи с гладким решением

2.4.1. Используемые обозначения.

2.4.2. Задача с доминирующей диффузией.

2.4.3. Задача с доминирующей конвекцией.

2.4.4. Задача с полным тензором диффузии.

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

Глава 3. Методы повышения эффективности вычислений

3.1. Параллелизация алгоритма.

3.1.1. Разбиение области по процессорам.

3.1.2. Параллелизация локальных шагов.

3.1.3. Параллелизация итерационного метода решения глобальной системы

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

3.2. Применение динамических сеток

3.2.1. Технология перестроения сетки

3.2.2. Переинтерполяция решения с одной сетки на другую

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

3.3. Двухуровневый переобуславливатель для задач диффузии с анизотропным неоднородным тензором диффузии.

3.3.1. Модельная задача, дискретизация и необходимые обозначения

3.3.2. Построение переобуславливателя.

3.3.3. Реализация алгоритма.

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

Глава 4. Решение прикладной трехмерной задачи о распространении ядерных загрязнений.

4.1. Постановка задачи.

4.2. Выбор расчетной тетраэдральной сетки.

4.3. Выбор расчетных схем

4.4. Результаты расчета гидравлического папора.

4.5. Результаты расчета распространения загрязнения.

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

Обычно модель долгосрочного распространения загрязнений предусматривает решение двух отдельных задач [20], [29]: нахождение фильтрационных потоков и моделирование процесса распространения. В настоящей работе фильтрационные потоки считаются стационариыми, бездивергентными и подчиняющимися закону Дарси. Распространение примеси в растворе описывается нестационарным уравнением конвекции-диффузии в пористой среде с дополнительными членами, отвечающими за химическое взаимодействие и радиоактивный распад. Основной целью диссертации является разработка, исследование и совершенствование численных методов решеиия именно второй задачи, но, в то же время, применяемые методы дискретизации диффузионного оператора пригодны и для задачи фильтрации.

Параметры реальных задач накладывают минимальные ограничения на коэффициенты диффузии и конвекции: тензоры проницаемости и диффузии могут быть полными, анизотропными, неоднородными (т.е. варьироваться при переходе от одной сеточной ячейки к другой), а поле конвективного переноса — сильно меняющимся в пространстве. Как следствие, в некоторых подобластях расчетной области доминирует диффузионный процесс, а в других — конвективный. Сложность геометрии расчетной области и наличие зон негладкости решения (например, вблизи источника или на границе раздела геологических слоев) предполагают использование неструктурированных локально сгущающихся, возможно, анизотропных расчетных сеток. Автором используются тетраэдральные сетки, как наиболее приспособленные для дискретизации сложных расчетных областей и локального измельчения. В сотрудничестве с Ю.В. Василевским и В.Н. Чугуновым создана программа на языке Fortran, позволяющая строить тетраэдральные сетки на основе гекса-эдральпых с последующим локальным измельчением.

Вышеизложенные свойства коэффициентов и расчетных сеток сужают класс возможных дискретизаций. Так, применение хорошо разработанных методов конечно-разностных дискретизаций [16, 19] для неструктурированных тетраэдральных сеток затруднительно. Наиболее перспективными аппроксимациями для рассматриваемого класса задач являются, па взгляд автора, методы конечных элементов (МКЭ)[15],[18] и методы конечных объемов (МКО)[35], в рамках которых будут рассмотрены некоторые подходы к построению вычислительных схем.

Разработка численных методов для уравнения конвекции-диффузии традиционно отталкивается от методов пространственной аппроксимации диффузионного оператора. Самым простым и популярным является МКЭ с пространством непрерывных кусочно-линейных функций [18], Pi-МКЭ. Недостатками этого метода являются немонотонность на произвольных неструктурированных сетках и отсутствие аппроксимации дискретных потоков на границах ячеек сетки. Несмотря на то, что возможно восстановление дискретных потоков на границах ячеек вспомогательной сетки (например, сетки Вороного), эти ячейки могут быть невыпуклы, содержать много граней и включать разрывы коэффициента диффузии. Стремление построить метод, приближающий диффузионные потоки на границах ячеек заданной сетки, привело к появлению других аппроксимаций. Метод смешанных конечных элементов (МСКЭ) [30] одновременно аппроксимирует как само решение, так и диффузионные потоки. Решение может иметь разрывы на границах элементов, диффузионные потоки находятся в пространстве Равьяра-Тома [56], благодаря чему обеспечивается непрерывность их нормальной составляющей на общих гранях двух элементов и, следовательно, локальная консервативность на ячейках сетки. К недостаткам МСКЭ в приложении к диффузионным задачам относится немонотонность схемы на неструктурированных и анизотропных сетках [28], [41]. Метод Галеркина с разрывными элементами (БС, 018соп{;тюи8 Са1ег1ап) [32], [25] является современной альтернативой вышеперечисленным методам в тех случаях, когда необходима гибкость при задании степени аппроксимирующего полинома для каждой ячейки сетки или используемые сетки являются неконформпыми. По своим свойствам метод близок к МСКЭ [27].

Метод конечных объемов (МКО) широко применяется в инженерных вычислениях. Существует множество вариаций МКО: методы для разных типов сеточных ячеек; методы со степенями свободы в вершинах сетки или в ячейках; методы с двухточечным или многоточечным шаблоном, аппроксимирующим поток через грань ячейки; методы с линейной или нелинейной аппроксимацией потока, и т.д. В применении к неструктурированным триан-гуляциям для МКО с двухточечной аппроксимацией диффузионного потока накладываются существенные ограничения на сетки (ортогональность сеток или сетки Делоне и Вороного [55],[37]), а для МКО с многоточечными шаблонами для потоков требуются сетки с регулярными ячейками [21], [34]. Кроме того, наличие полного тензора диффузии существенно усложняет задачу построения линейной схемы МКО. В диссертационной работе предложен новый монотонный нелинейный МКО для произвольных сеток и тензоров диффузии (или проницаемости для фильтрационных задач) [11]. Метод использует идеи, опубликованные в работах [52],[53] для двумерных задач диффузии. Под монотонностью метода здесь и далее подразумевается свойство неотрицательности получаемого численного решения при соответствующих граничных условиях и правой части, обеспечиваемое монотонностью матрицы ([3], стр.269) аппроксимации диффузионного оператора, аналогично работам [52], [53]. В случае однородных краевых условий типа Дирихле такая трактовка понятия монотонности соответствует принципу максимума с диагональным преобладанием по столбцам в матрице аппроксимации, сформулированному А.А.Самарским и П.Н.Вабищевичем в [16](стр.63). Неотрицательность численного решения приобретает исключительную важность в задачах, где требуется учет химических взаимодействий, что отмечено в работе [36]. Наличие отрицательной концентрации в расчетной области может приводить к неверному расчету химических реакций и неконсервативности модели.

Переход от задачи диффузии к задаче конвекции-диффузии сопряжен с дополнительными трудностями при построении устойчивых схем в случае доминирующей конвекции. Во-первых, для стабилизации схем могут применяться методы регуляризации [16]. Классическим примером является Pi-МКЭ в сочетании с алгоритмом SUPG (Streamline Upwinding Petrov Galerkin) [42]. В случае МСКЭ и МКО используется противопотоковая аппроксимация конвективного члена. Во-вторых, для аппроксимации конвективного оператора возможно использование схем сквозного счета, основанных на методах Годунова [4], Лакса и Фридрихса [49] и др. Примерами таких схем являются методы, использующие разрывные конечные элементы (РКЭ) [58], [33], а также метод коррекции потока [44], обеспечивающий монотонность и минимальную численную диффузию, который, однако, не применим к случаю полных и неоднородных тензоров диффузии. Эти методы позволяют аппроксимировать негладкие решения с малой численной диффузией и без осцилляций.

Для аппроксимации по времени нестационарных задач конвекции-диффузии могут использоваться явные схемы, неявные схемы и методы расщепления [14]. При этом явные и неявные схемы используют одну и ту же пространственную аппроксимацию конвективных и диффузионных членов (МКЭ, МСКЭ, МКО; сравнение методов на двумерной тестовой задаче проведено в [28]), а методы расщепления допускают построение разных пространственных аппроксимаций отдельно для каждого из физических операторов. Так, например, МСКЭ может использоваться для аппроксимации диффузионных потоков, РКЭ - для аппроксимации решения [58]. Доусон [33], Аккерер и др. [23] предлагают схемы расщепления с явной аппроксимацией конвекции и неявной — диффузии. Вохралик [61] использует МКО для конвективного члена и иеконформные КЭ или МСКЭ для диффузионного оператора в сочетании с полностью неявной схемой. Схемы расщепления позволяют для каждого из пространственных операторов выбирать наиболее адекватный метод аппроксимации, поэтому именно они рассматриваются в работе как наиболее приспособленные для моделирования процессов переноса примесей в пористых средах.

Рассматриваемый класс трехмерных проблем относится к разряду больших задач, требующих больших вычислительных мощностей, машинной памяти, расчеты занимают дни и недели. Для повышения эффективности вычислений применяются различные технологии: в работах [33], [23] предлагаемые схемы расщепления допускают применение разных шагов по времени для конвекции и диффузии, а также локальных шагов по времени в подобластях; В.И.Агошковым и др. [24] предложен метод декомпозиции области для процессов конвекции-диффузии в неоднородных средах, который позволяет разбить задачу на несколько подзадач в подобластях, решая только задачу конвекции в подобластях с сильным преобладанием коивекции над диффузией. В настоящей работе предложен ряд методов ускорения вычислений. Па-раллелизация алгоритма [2],[10] дает возможность решать задачи на сетках, содержащих несколько миллионов тетраэдров. Технология автоматического динамического перестроения (локального сгущения и разгрубления) сеток сокращает время работы программы благодаря меньшему числу неизвестных и обеспечивает высокую точность в областях быстрого изменения решения [9]. Двухуровневый переобуславливатель для задач диффузии в пластовых средах [46],[47], разработанный и исследованный в сотрудничестве с Ю.А. Кузнецовым, убирает зависимость скорости сходимости метода сопряженных градиентов от коэффициентов диффузии и шага сетки по вертикальной оси.

Работы с результатами численного моделирования реальных объектов захоронения загрязнений редко встречаются в печати. В тематическом выпуске журнала Computational Geosciences (т.8, №2, 2004) опубликованы результаты расчетов задания Couplex, предложенного французским агентством по ядерным отходам. Задание состояло из трех задач: трехмерного моделирования мелкомасштабных процессов конвекции-диффузии с учетом химических реакций в малой окрестности захоронения ядерных отходов; двумерного моделирования крупномасштабного переноса загрязнения от источника; совмещения мелкомасштабной модели источника загрязнений с крупномасштабной моделью двумерного переноса. П.И.Чаловым и др. в статье [20] представлены результаты расчета переноса загрязнений от хвостохрапилища карабал-тинского гидрометаллургического завода по переработке урапосодержащих руд (Кыргызстан). В этой задаче поле конвективных потоков было двумерным, для дискретизации использовались гексаэдральпые сетки. В настоящей работе проведено численное моделирование новой полностью трехмерной задачи крупномасштабного переноса загрязнения па неструктурированных тетраэдральных сетках. Задание представляет большой интерес с методической точки зрения, т.к. расчет затрудняют следующие факторы: анизотропия (2 порядка) и неоднородность (до 10 порядков) тензора проницаемости; сильная анизотропия области, и, как следствие, расчетной сетки (диаметр области — 55км, глубина — 500м); неоднородность и анизотропия полного тензора диффузии; выклинивания геологических пластов; необходимость прогноза на большой период времени (сотни тысяч лет).

Основные результаты диссертационной работы:

1. Разработано новое семейство монотонных методов конечных объемов для аппроксимации диффузионных задач с полным и неоднородным тензором диффузии на неструктурированных тетраэдральных сетках.

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

3. Предложены и протестированы новые алгоритмы повышения эффективности вычислений: решение задач переноса на динамически перестраиваемых расчетных сетках; двухуровневый переобуславливатель для диффузионных задач со скачками анизотропного тензора диффузии.

4. Проведено численное моделирование трехмерной задачи с реальными параметрами.

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

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

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

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

Автор выражает глубокую признательность научному руководителю Ю.В.Василевскому за постановку задачи и полезные обсуждения; Ю.А.Кузнецову за предоставленную возможность участвовать в разработке и исследовании нового переобуславливателя; Д.А. Святскому, В.Н.Чугунову, и К.Н.Липникову за плодотворные дискуссии по возникавшим в ходе написания диссертации вопросам; коллективу ИВМ РАН за создание благоприятной рабочей атмосферы. Такэюе автор благодарен первому научному руководителю И.В.Новоэюилову, пробудившему интерес к научной работе.

Выводы по главе 4

В главе проведено численное моделирование прикладной задачи Andra 3D, включающее построение тетраэдральной сетки и раздельное решение задач фильтрации и конвекции-диффузии. Рассмотрены три подхода к построению тетраэдральных сеток, выбран наиболее адекватный метод для данной задачи и создана сетка, содержащая 2.5 миллиона ячеек. Найдено поле фильтрационных потоков и гидравлический напор в условиях сильной анизотропии и неоднородности тензора фильтрации и расчетной сетки.

Решение задачи переноса также осложнялось сжатием сетки по вертикали, анизотропией и неоднородностью тензора диффузии, разными временными масштабами процессов в верхних и нижних геологических пластах. В результате расчетов получен прогноз времени и места выхода на земную поверхность загрязненных вод с концентрациями радиоактивного Йода 129, равными Ю-14, Ю-12, Ю~10кг/м3 . Показана динамика распространения загрязнения в трехмерной области с течением времени, вплоть до 350000 лет.

Заключение

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

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

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

1. Васпиев К. С., Дмитриев Н. М., Каневская Р. Д., Максимов В. М. Подземная гидромеханика. — М.-Ижевск: Институт компьютерных исследований, 2005. - 496 с.

2. Василевский Ю. В., Капырин И. В. Параллельное трехмерное моделирование распространения примесей в пористых средах // Матричные методы и технологии решения больших задач. — М.: ИВМ РАН, 2005. — С. 33-50.

3. Воеводин В. В., Кузнецов Ю. А. Матрицы и вычисления. — М.: Наука, 1984.- 320 с.

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

5. Ентов В. М. Теория фильтрации // Соросовский образовательный оюур-нал. 1998. - № 2. - С. 121-128.

6. Интернет-ресурс: http://sourceforge.net/projects/ani2d/.

7. Интернет-ресурс: http://sourceforge.net/projects/ani3d/.

8. Каневская Р. Д. Математическое моделирование гидродинамических процессов разработки месторождений углеводородов. — М.-Ижевск: Институт компьютерных исследований, 2003. — 128 с.

9. Капырин И. В. Об использовании динамических сеток при решении задач конвекции-диффузии // Труды математического центра им. Н.И. Лобачевского. 2006. - Т. 31. - С. 166-175.

10. Капырин И. В. Семейство монотонных методов численного решения трехмерных задач диффузии на неструктурированных тетраэдральных сетках // Доклады Академии Наук. 2007. — Т. 614, № 5. — С. 588-593.

11. Кузнецов Ю. А. Итерационные методы в подпространствах. — М.: ОВМ АН СССР, 1984. 136 с.

12. Лебедев В. ИАгошков В. И. Операторы Пуанкаре-Стеклова и их приложения в анализе. М.: ОВМ АН СССР, 1983.- 184 с.

13. Марчук Г. И. Методы расщепления. — М.: Наука, 1988. — 264 с.

14. Марчук Г. И., Агошков В. И. Введение в проекционио-сеточные методы. М.: Наука, 1981. - 416 с.

15. Самарский А. А., Вабищевич П. Н. Численные методы решения задач конвекции-диффузии. — М.: Эдиториал УРСС, 1999.— 248 с.

16. Седов Л. И. Механика сплошных сред. — М.: Наука, 1970. — Т. 1. — 492 с.

17. Сьярле Ф. Метод конечных элементов для эллиптических задач. — М.: Мир, 1980.- 511 с.

18. Толстых А. И. Компактные разностные схемы и их применение в задачах аэрогидродинамики. — М.: Наука, 1990.— 230 с.

19. Aavatsmark I., Barkve Т., Вое О., Mannseth Т. Discretization on unstructured grids for inhomogencous, anizotropic media, part I: Derivation of the methods // SIAM. J. Sci. Comput.- 1998.- Vol. 19, no. 5.-Pp. 1700-1716.

20. Achdou Y., Jaffre J., Svyatski D., Vassilevski Y. Numerical simulation of unsteady 3d flows on anisotropic grids // Transactions of French-Russian Lia-pounov Institute for Applied Mathematics and Computer Science. — Moscow, 2003.-Vol. 4,- Pp. 107-124.

21. Ackerer P., Younes A., Mose R. Modeling variable density flow and solute transport in porous medium: 1. Numerical model and verification // Transport in Porous Media. 1999. - Vol. 35. - Pp. 345-373.

22. Agoshkov V., Gervasio P., Quarteroni A. Optimal control in heterogeneous domain decomposition methods for advection-diffusion equations // Mediterranean Journal of Mathematics. — 2006. — Vol. 3. — Pp. 147-176.

23. Arnold В., Brezzi F., Cockburn В., Marini L. Unified analysis of discontinuous Galerkin methods for elliptic priblems // SIAM J.Numer.AnaL — 2002.- Vol. 39, no. 5,- Pp. 1749-1779.

24. Bansch E. Local mesh refinement in 2 and 3 dimensions // IMPACT of Computing in Science and Engrg. — 1991. — Vol. 3. — Pp. 181-191.

25. Bastian P., Lang S. Couplex benchmark computations obtained with the software toolbox UG // Computat. Geosci 2004. - no. 8.- Pp. 125-147.

26. Bernard-Michel G., Le Potier C., Beccantini A. et al The Andra Couplex 1 test case: Comparisons between finite element, mixed hybrid finite element and finite volume discretizations // Computat. Geosci.— 2004.— Vol. 8, no. 2.- Pp. 187-201.

27. Bourgeat A., Kern M., Schumacher S., Talandier J. The COUPLEX test cases: Nuclear waste disposal simulation // Computat. Geosci— 2004.— Vol. 8, no. 2. Pp. 83-98.

28. Brezzi F., Fortin. M. Mixed and hybrid finite element methods. — New York: Springer Verlag, 1991.- 350 p.

29. Chavent G., Jaffre J. Mathematial models and finite elements for reservoir simulation. — North Holland, Amsterdam, 1986. — 388 p.

30. Cockburn B., Karniadakis G. E., Shu C. Discontinuous Galerkin methods, theory, computation and applications // Lecture Notes in Computational Science and Engineering. — 2000.— Vol. 11.

31. Dawson C. High resolution upwind-mixed finite element methods for advec-tion-diffusion equations with variable time-stepping // Numerical methods for Partial Differential equations. — 1995. — Vol. 11, no. 5.— Pp. 525-538.

32. Droniou J., Eymard R. A mixed finite volume scheme for anisotropic diffusion problems on any grid // arXiv.org: math/0604010. — 2006.

33. Eymard R., Gallouet T., Herbin R. Finite volume methods // Handbook of numerical analysis / Ed. by P. Ciarlet, J. Lions.— Amsterdam, North Holland, 2000. Vol. 7. - Pp. 713-1020.

34. Herbin R., Labergerie 0. Finite volume schemes for elliptic and elliptic-hyperbolic problems on triangular meshes // Comput. Methods Appl. Mech. Engrg. 1997. - no. 147. - Pp. 85-103.

35. Hoteit H. Simulation d'écoulement et de transports de polluants en milieu poreux: Application à la modélisation de la sûreté des dépôts de déchets radioactifs: Ph.D. thesis / l'Université de Rennes 1.— 2002.

36. Hoteit H., Ackerer P., Mose R. Nuclear waste disposal simulations: Couplex test cases. // Comput. Geosci. 2004. - Vol. 8, no. 2.- Pp. 99-124.

37. Hoteit H., Ackerer P., Mose R. et al. New two-dimensional slope limiters for discontinuous Galerkin methods on arbitrary meshes // Int. J. Numer. Meth. Engng.-20Qi.-Wol 61.- Pp. 2566-2593.

38. Hoteit H., Mose R., Philippe B. et al. The maximum principle violations of the mixed-hybrid finite-element method applied to diffusion equations // Numer. Meth. Engng. 2002. - Vol. 55, no. 12,- Pp. 1373-1390.

39. Hughes T. J. R., Brooks A. N. A multidimensional upwind scheme with no crosswind diffusion // Finite Element Methods for Convection Dominated Flows. 1979. - Vol. 34. - Pp. 19-35.

40. Kaporin I. High quality preconditioning of a general symmetric positive definite matrix based on its UTU+UTR+RTU-dQcomposit\on // Numer.Linear Algebra Appl. 1998. - Vol. 5. - Pp. 483-509.

41. Kuzmin D., Turek S. Flux correction tools for finite elements // J. Comput. Phys. 2002. - Vol. 175, no. 2. - Pp. 525-558.

42. Kuznetsov Y. A. Two-level preconditioners with projectors for unstructured grids // Russian Journal of Numerical Analysis and Mathematical Modelling. 2000. - Vol. 15, no. 3-4. - Pp. 247-255.

43. Kuznetsov Y. A. Domain decomposition preconditioner for anisotropic diffusion // Lecture Notes in Computational Science and Engineering: Domain Decomposition Methods in Science and Engineering XVII. — Vol. 60. — Springer, 2007.

44. Kuznetsov Y., Lipnikov K. An efficient iterative solver for a simplified poroe-lasticity problem // East-West Journal of Numerical Mathematics. — 2000. — Vol. 8, no. 3.-Pp. 207-221.

45. Lax P. D. Shock waves and entropy // Contributions to nonlinear functional analisys. — New York: Academic Press, 1971. — Pp. 603-634.

46. Leij F. J., Dane J. H. Analitical solutions of the one-dimensional advec-tion equation and two- or three-dimensional dispersion equation // Water Resources Research.- 1990,- Vol. 26, no. 7.- Pp. 1475-1482.

47. Leij F. J., Dane J. H., van Genuchten M. T. Mathematical analysis of one-dimensional solute transport in a layered soil profile // Soil. Sci. Soc. Am. — 1991.-no. 55.-Pp. 944-953.

48. Le Potier C. Schema volumes finis monotone pour des operateurs de diffusion fortement anisotropes sur des maillages de triangle non structures // C.R.Acad. Sei. Paris, Ser. 1341. 2005. - Pp. 787-792.

49. Lipnikov K., Vassilevski Y. Parallel adaptive solution of 3d boundary value problems by Hessian recovery // Comp.Methods Appl.Mech.Engnr.— 2003.-Vol. 192.-Pp. 1495-1513.

50. Mishev I. D. Finite volume methods on Voronoi meshes // Numerical Methods for Partial Differential Equations.— 1998.— Vol. 12, no. 2.— Pp. 193-212.

51. Raviart P., Thomas J. A mixed hybrid finite element method for the second order elliptic problem // Lectures Notes in Mathematics.— 1977.— Vol. 606.-Pp. 292-315.

52. Saad Y. Iterative methods for sparse linear systems, Second Edition.— Philadelphia, PA: SIAM, 2000. 460 p.

53. Siegel P., Mose R., Ackerer P., Jaffre J. Solution of the advection-diffusion equation using a combination of discontinuous and mixed finite elements // Internat. J. Numer. Methods Fluids. 1997. - Vol. 24. - Pp. 595-613.

54. Stuben K. Algebraic multigrid: experience and comparisons // Applied Math, and Comp. 1983. - Vol. 13, no. 3-4. - Pp. 419-451.

55. Vassilevski Y. A hybrid domain decomposition method based on aggregation // Numer. Linear Algebra AppL— 2004. — no. 11.— Pp. 327-341.

56. Vohralik M. Numerical Methods for nonlinear elliptic and parabolic equations: Ph.D. thesis / University Paris XI, Chech Technical University of Prague. 2004.

57. Williams R. Performance of dynamic load balancing algorithms for unstructured mesh calculations. // Concurrency: Practice and Experience. — 1992. — no. 3. Pp. 457-481.