Содержание к диссертации
Введение
Глава 1 Уравнения фильтрации двухфазной жидкости 17
1.1 Математическая постановка задачи 17
1.1.1 Основные физические параметры 17
1.1.2 Уравнения неразрывности и закон Дарси 18
1.1.3 Капиллярное давление 20
1.1.4 Система, как уравнения первого порядка 20
1.1.5 Краевые условия 22
1.2 Смешанная обобщенная формулировка 26
Глава 2 Численное решение задачи 30
2.1 Пространственная аппроксимация 31
2.1.1 Элементы RT[0], двумерная задача 31
2.1.2 Элементы RT[0], трехмерная задача 38
2.1.3 Итоговая сеточная задача 43
2.2 Интегрирование по времени 44
2.2.1 Явная схема предиктор-корректор 44
2.2.2 Алгоритм решения задачи фильтрации 45
2.2.3 Выполнение интегрального баланса 47
2.3 Решение седловой задачи 48
2.3.1 Свойства СЛАУ 48
2.3.2 Переобусловленный метод сопряженных градиентов 49
2.4 Построение переобуславливателя 50
2.4.1 Решение одномерной спектральной задачи 51
2.4.2 Двумерный случай 53
2.4.3 Трехмерный случай 55
2.4.4 Сравнение со стандартным оператором Лапласа . 58
Глава 3 Моделирование процесса фильтрации 60
3.1 Моделирование скважин 61
3.2 Двумерная задача 62
3.2.1 Сходимость 62
3.2.2 Плановая задача 63
3.2.3 Вертикальная задача 69
3.3 Трехмерная задача 74
3.3.1 Сходимость 74
3.3.2 Примеры скважин 75
3.4 Влияние капиллярного давления 84
Глава 4 Параллелизация алгоритма 86
4.1 Двумерный случай 87
4.1.1 Распределение данных по процессам 87
4.1.2 Реализация алгоритма 88
4.1.3 Результаты ускорения вычислений 91
4.1.4 Архитектуры с распределенной и общей памятью . 95
4.2 Трехмерный случай 98
4.2.1 Распределение данных по процессам 98
4.2.2 Реализация алгоритма 99
4.2.3 Результаты ускорения вычислений 102
4.2.4 Сетевой закон Амдала 103
Заключение 107
Литература 109
- Уравнения неразрывности и закон Дарси
- Элементы RT[0], трехмерная задача
- Вертикальная задача
- Архитектуры с распределенной и общей памятью
Введение к работе
Выбор темы исследования обусловлен тем, что в современных условиях совместная фильтрация несмешивающихся жидкостей является одним из ключевых разделов в подземной гидродинамике, а определяющим инструментом по получению новых знаний выступает математическое моделирование.
Актуальность работы определяется необходимостью создания алгоритмов и средств, обеспечивающих эффективное решение задачи совместной фильтрации несмешивающихся жидкостей методами математического моделирования на вычислительных системах терафлопного диапазона.
В результате активного развития вычислительной техники производительность существующих на данный момент многопроцессорных систем измеряется уже сотнями терафлопс. Однако новая архитектура предъявляет особые требования к алгоритмам и их реализациям. Большое количество процессов требует высокой степени масштабируемости алгоритма, многоядерные узлы -особенного внимания к оптимизации доступа к памяти. Чтобы эти расчеты были эффективны, на таких системах требуется расширенная параллельная реализация, соответствующая кластерной структуре вычислительной системы с многопроцессорными узлами.
Целью работы является построение «технологической цепочки», соответствующей процессам совместной фильтрации несмешивающихся жидкостей, и создание алгоритмов и компонентов вычислительной среды, обеспечивающих высокую масштабируемость и возможность эффективного решения поставленной задачи на многопроцессорных вычислительных системах терафлопного диапазона. Схематически эта цепочка выглядит следующим образом: от изучаемого явления - к его математической модели; затем - к аппроксимационной модели; далее - к эффективному численному алгоритму; программе, реализующей этот алгоритм; и, наконец, к вычислениям на суперкомпьютерах, анализу полученных результатов, уточнению в случае необходимости математической модели.
Научная новизна.
Реализована и экспериментально обоснована численная модель Маскета-Леверетта в смешанной обобщенной постановке в терминах «скорость-давление-насыщенность» с произвольным числом скважин, моделирование которых не требует сгущения сетки в прискважинных зонах и соразмерности диаметра скважин шагу сетки.
Создан эффективный алгоритм для обращения сеточного седлового оператора, возникающего в смешанной постановке и являющегося аналогом оператора Лапласа, используемый в качестве переобуславливателя в итерационном методе сопряженных градиентов при решении эллиптического уравнения второго порядка с переменными коэффициентами.
- Разработан и изучен алгоритм распараллеливания двумерной и трехмер
ной задач совместной фильтрации на суперкомпьютерах, удовлетворяющий
кластерной структуре с использованием мрі технологии в рамках предложен
ной модели.
Практическая значимость диссертационной работы заключается в возможности использования разработанного комплекса программ для интенсификации разработки месторождений, а также в универсальности предлагаемой методики распараллеливания и ее применимости к другим задачам механики сплошной среды.
Достоверность и обоснованность полученных результатов подтверждается соответствием рассматриваемой модели фундаментальным законам сохранения. Для предложенной численной модели показано, что используемая аппроксимация в комбинации с явной схемой обеспечивают выполнение уравнения баланса; приводятся расчеты, демонстрирующие её аппроксимационные свойства и адекватность известным явлениям, сопутствующим процессу вытеснения нефти водой. Эти результаты хорошо согласуются с работами А.Н. Коновалова, Б.Н. Четверушкина, В.И. Дробышевича. Для разработанной методики распараллеливания приводятся графики, подтверждающие её эффективность.
Основные результаты, выносимые на защиту:
Реализация и экспериментальное обоснование численной модели Маскета-Леверетта в смешанной обобщенной постановке в терминах «скорость-давление-насыщенность» с произвольным числом скважин.
Эффективный алгоритм для обращения сеточного седлового оператора, возникающего в смешанной постановке и являющегося аналогом оператора Лапласа, используемый в качестве переобуславливателя в итерационном методе сопряженных градиентов при решении эллиптического уравнения второго порядка с переменными коэффициентами, которое является определяющим в задаче фильтрации.
Параллельный алгоритм и программный комплекс для решения двумерной и трехмерной задач совместной фильтрации на суперкомпьютерах, удовлетворяющий кластерной структуре и использующий мрі технологии в рамках предложенной модели.
Апробация диссертации. Результаты работ по диссертации докладывались на научно-исследовательских семинарах ИВМиМГ СО РАН, на Конференциях молодых ученых, на Всероссийских конференциях по вычислительной математике, на XLV Международной научной студенческой конференции, на Международной конференции «Дифференциальные уравнения. Функциональные пространства. Теория приближений».
Работа поддержана проектом Президиума СО РАН и грантом РФФИ.
Публикации. По теме диссертации опубликовано 9 печатных работ, из них 3 статьи в ведущих рецензируемых журналах и изданиях, перечень которых
определен Высшей аттестационной комиссией.
Личный вклад автора. Все выносимые на защиту результаты принадлежат лично автору. Автор принимал активное участие в постановке двумерной и трехмерной задач совместной фильтрации. Реализовал и обосновал численную модель, разработал параллельный алгоритм, который удовлетворяет кластерной структуре и использует MPI технологи. Провел эксперименты, интерпретировал их результаты, показал высокую масштабируемость программного комплекса. Автор разработал алгоритм для обращения сеточного седлового оператора, возникающего в смешанной постановке и являющегося аналогом оператора Лапласа, используемый в качестве переобуславливателя в итерационном методе сопряженных градиентов при решении эллиптического уравнения второго порядка с переменными коэффициентами.
Структура и объем работы. Диссертация состоит из введения, четырех глав, заключения и списка литературы. Работа содержит 33 рисунка, 10 таблиц, 71 наименование библиографии. Полный объем диссертации составляет 118 страниц.
Уравнения неразрывности и закон Дарси
Достоверность и обоснованность полученных результатов подтверждается соответствием рассматриваемой модели фундаментальным законам сохранения. Для предложенной численной модели показано, что используемая аппроксимация в комбинации с явной схемой обеспечивают выполнение уравнения баланса; приводятся расчеты, демонстрирующие её аппроксимационные свойства и адекватность хорошо известным явлениям, сопутствующим процессу вытеснения нефти водой. Эти результаты хорошо согласуются с работами Дробышевича, Коновалова, Четве-рушкина. Для разработанной методики распараллеливания приводятся графики, подтверждающие её эффективность.
Перейдем к краткому описанию содержания диссертации. Диссертация состоит из введения, четырех глав, заключения и списка цитируемой литературы. Для удобства чтения каждая глава предваряется кратким введением. Заключение содержит резюме о полученных результатах. Ссылки на первоисточники даны во введении. В основной части текста упоминаются лишь работы, содержащие некоторые конкретные факты, используемые для доказательств утверждений. Каждая глава разделена на пункты с трехиндексными номерами. В диссертации принята сквозная трехипдексная нумерация формул, теорем, лемм и ссылок на них. Первый индекс соответствует номеру главы, второй - номеру пункта главы, третий - номеру формулы или утверждения данной главы. Работа содержит 33 рисунка, 10 таблиц, 71 наименование библиографии. Полный объем диссертации составляет 118 страниц.
Во введении сформулирована цель работы, обоснованы актуальность темы и научная новизна, кратко описано содержание работы.
В первой главе формулируется математическая постановка модели Маскета-Леверетта для задачи двухфазной фильтрации несжимаемой жидкости, в терминах «скорость-давление-насыщенность». Приводятся основные соотношения, определяющие процесс фильтрации, формулируется задача о вытеснении нефти водой в виде системы уравнений первого порядка с соответствующими краевыми условиями. Особое внимание уделено обсуждению условий непротекания фаз при наличии силы тяжести. От системы дифференциальных уравнений осуществляется переход к слабой постановке в интегральных уравнениях.
Вторая глава посвящена аппроксимациям интегральных тождеств в слабой постановке по смешанному методу конечных элементов для двумерного и трехмерного случаев. При этом векторные поля аппроксимируются элементами Равьяра-Тома наименьшей степени на прямоугольной сетке, а скалярные функции (давление и насыщенность) ищутся в пространстве кусочно-постоянных функций. Предлагается алгоритм, основанный на решении сеточных уравнений для скорости и давления на одном шаге интегрирования по времени. Указанные уравнения представляют собой линейную алгебраическую систему седлового типа, решение которой осуществляется методом сопряженных градиентов для дополнения Шура с переобусловливателем, допускающим разделение переменных и вычисления по формулам быстрого дискретного преобразования Фурье (см. [26, 61]). Отметим, что система для дополнения Шура соответствует уравнению для давления после исключения суммарной скорости и является аналогом оператора Лапласа. В качестве разностной схемы по времени используется явная схема типа предиктор-корректор второго порядка аппроксимации, но только с одним вычислением правой части на шаге интегрирования (см. [8]).
В третье главе приводятся результаты численных экспериментов, демонстрирующие аппроксимационные свойства численной модели и её адекватность хорошо известным явлениям, сопутствующим процессу вытеснения нефти водой. Также в этой главе приводится способ моделирования нагнетательных и добывающих скважин. Особое внимание уделено исследованию поведения скважин и групп скважин, определению остаточных запасов, застойных зон на конкретные моменты времени, оценке запасов по пластам и в целом по залежи, выбору оптимальных режимов работы скважин, а также сокращения общего числа скважин. Все это позволяет при освоении нефтяного месторождения проводить анализ рисков и своевременно минимизировать их, совершенствовать технологии, обосновывать стратегии доразработки, улучшать показатели добычи.
В четвертой главе для эффективного решения двумерной и трехмерной задач фильтрации на суперкомпьютерах предлагается алгоритм распараллеливания, удовлетворяющий кластерной структуре с использованием MPI [65, 66, 67, 68] технологии, реализованный на языке Fortran. Данный алгоритм равномерно распределяет все данные по процессам и одинаково загружает их в каждый момент времени выполнения программы. Представлены результаты экспериментов, которые проводились в Сибирском суперкомпьютерном центре (ССКЦ) [70] на кластере нкс-зот. Отметим, что для передачи данных между узлами используется прямой доступ к памяти с помощью RDMA [69] технологии через сетевые карты InfiniBand [64], которые имеют низкую задержку и высокую пропускную способность. Приведены зависимости производительности от числа процессов при различных размерностях задачи, на основании которых можно сделать вывод об оптимальном числе процессов. Также в главе сравниваются архитектуры с распределенной и общей памятью. Показывается, что выгодно строить параллельные системы с использованием обоих видов архитектур. Приводится закон Амдала и его расширенный вариант, учитывающий потери времени на межпроцессорный обмен сообщениями. Ведь даже если задача обладает идеальным параллелизмом, то узким местом (bottleneck) является железо, а именно сетевые карты. В заключении данной диссертационной работы дано краткое изложение основных результатов.
Элементы RT[0], трехмерная задача
Каждая из этих матриц состоит из трех блоков, первый хранится в направлении х, второй - в направлении у, а третий - в направлении z. Для того, чтобы эффективно реализовать умножение вектора на матрицу, соответствующую одному из направлений, мы должны распределить вектор по процессам в этом же направлении (перенумеровать его). Изначально мы храним скалярные функции в направлении х, а каждую компоненту (ж, у, z) векторных функции храним в соответствующем направлении. Для того, чтобы перераспределить скалярные функции по процессам в направлениях у и z, используем стандартную процедуру MPI_ALLTOALL. Результатом умножения на матрицу В является вектор, структура которого соответствует структуре векторного поля v. Матрица D(s) состоит из трех независимых блоков: причем каждый из процессов содержит независимые трехдиагональные блоки матриц Dx(s), Dy{s) и Dz(s). Тогда действие матриц .JD 1(s), D l(s) и D l(s) на, соответственно, х, у и z компоненты вектора реализуется прогонками на каждом из процессов. Умножение на матрицу
В осуществляется аналогично умножению на матрицу В, только при этом процессы обмениваются данными в обратном порядке по сравнению с умножением на В.
Рассмотрим вопрос о параллельной реализации обращения на векторе невязки переобуславливателя (2.4.10). С точки зрения работы с памятью указанная процедура аналогична описанному выше умножению на мат-рицы В, D (s) и В . Сначала в соответствии с первым разложением находятся коэффициенты Фурье вектора / в направлениях х и у. Отметим, что, как и в случае с матрицами для эффективной реализации метода Фурье по одному из направлений, мы должны распределить вектор по процессам в этом же направлении (перенумеровать его). Затем полученный вектор перераспределяется по процессам в направлении z. Следующий шаг состоит в решении трехдиагональных систем (2.4.16). И, наконец, в соответствии со вторым разложением (2.4.11), перераспределяя вектор сначала в направлению у, а затем в направлении ж, производим обратное преобразование Фурье.
Число процессов в нашей постановке не должно превосходить минимальную размерность задачи. Обычно для задач фильтрации это размерность по оси z. Это связано с тем, что в случае неоднородности области, требуется использовать «слоистый» переобуславливатель. Данный переобуславливатель отличается от приведенного в пункте 2.4.3 тем, что не требует по оси z однородности данных. Чтобы применить данный вид переобуславливателя, требуется осреднить область по каждому сечению XY.
Отметим, что предложенный в работе переобуславливатель легко обобщается для случая «слоистой» области, поскольку преобразования Фурье требуются только для осей х и /, а для оси z используется прогонка. Вместо конкретной матрицы Az, описанной в пункте 2.4.3, будет применяться произвольная матрица, но матрица эта все равно останется трехдиагональной.
Ниже приводятся расчеты коэффициентов ускорений и загруженности (рисунки 4.9, 4.10 и таблицы 4.4, 4.5)) на последовательностях сеток 32 х 32 х 4, 64 х 64 х 8, 128 х 128 х 16, 256 х 256 х 32 и 512 х 512 х 64 с использованием до 64 процессов кластера НКС-30.
Из приведенных результатов видно, что для самой густой сетки 512 х 512 х 64 имеется ухудшение результатов масштабируемости для максимального числа процессов. Вероятнее всего это связано с тем, что данный суперкомпьютер, на котором производились эксперименты имеет пропускную способность не более 20Gbit/s (Double Data Rate). В то время как, например, суперкомпьютер «Ломоносова» [71] (число вычислительных узлов - 5130, число процессоров - 10260, число процессорных ядер - 44000), занимавший 17 место в ТОР500 на 2010 год, имеет InfiniBand карты с пропускной способностью 40Gbit/s (Quad Data Rate), т.е. пересылает данные в два раза быстрее, что существенно сказывается на масштабируемости при больших объемах данных. Для лучшего понимания данного факта, рассмотрим сначала общеизвестный закон Амдала, а далее его модификацию - сетевой закон Амдала.
Одной из главных характеристик параллельных систем является ускорение Sp, которое уже было определено в пункте 4.1.3. Пусть W — Ws + WP7 где W - общее число операций в задаче, Wp - число операций, которые можно выполнять параллельно, a Ws — число последовательных (нераспараллеливаемых) операций. Обозначим также через t время выполнения одной операции. Тогда получаем известный закон Амдала ([44]), который показывает, что ускорение, которое может быть получено на вычислительной системе из р процессов, по сравнению с последовательным решением не будет превышать величины
Вертикальная задача
Из приведенных результатов видно, что для самой густой сетки 512 х 512 х 64 имеется ухудшение результатов масштабируемости для максимального числа процессов. Вероятнее всего это связано с тем, что данный суперкомпьютер, на котором производились эксперименты имеет пропускную способность не более 20Gbit/s (Double Data Rate). В то время как, например, суперкомпьютер «Ломоносова» [71] (число вычислительных узлов - 5130, число процессоров - 10260, число процессорных ядер - 44000), занимавший 17 место в ТОР500 на 2010 год, имеет InfiniBand карты с пропускной способностью 40Gbit/s (Quad Data Rate), т.е. пересылает данные в два раза быстрее, что существенно сказывается на масштабируемости при больших объемах данных. Для лучшего понимания данного факта, рассмотрим сначала общеизвестный закон Амдала, а далее его модификацию - сетевой закон Амдала.
Одной из главных характеристик параллельных систем является ускорение Sp, которое уже было определено в пункте 4.1.3. Пусть W — Ws + WP7 где W - общее число операций в задаче, Wp - число операций, которые можно выполнять параллельно, a Ws — число последовательных (нераспараллеливаемых) операций. Обозначим также через t время выполнения одной операции. Тогда получаем известный закон Амдала ([44]), который показывает, что ускорение, которое может быть получено на вычислительной системе из р процессов, по сравнению с последовательным решением не будет превышать величины свойствами алгоритма, а гт — техническую составляющую, которая зависит от соотношения быстродействия процессора и аппаратуры сети. Таким образом, для повышения скорости вычислений следует воздействовать на обе составляющие коэффициента деградации. Эти коэффициенты определяются множеством факторов: алгоритмом задачи, размером данных, реализацией функций обмена библиотеки MPI, использованием разделяемой памяти и, конечно, техническими характеристиками коммуникационных сред и их протоколов (именно поэтому на суперкомпьютерах используются сетевые карты InfiniBand, которые имеют низкую задержку и высокую пропускную способность).
Причем, если г — 0, то Sp — р, т.е. стремится к максимуму (см. 4.2.2). Следовательно, сетевой закон Амдала должен быть основой оптимальной разработки алгоритма и программирования задач, предназначенных для решения на многопроцессорных вычислительных системах. Перечислим основные результаты, полученные в диссертации. 1. На основе смешанного метода конечных элементов построена численная модель фильтрации двухфазной несжимаемой жидкости в терминах «скорость-давление-насыщенность». Данный подход позволил избежать в модели Маскета-Леверетта проблему вырождаемости краевых условий для насыщенности. 2. Приведен способ моделирования скважин как в двумерном, так и в трехмерном случае, не требующий сгущения сетки в прискважин-ных зонах и соразмерности диаметра скважин шагу сетки. 3. Создан эффективный алгоритм для обращения сеточного седлово-го оператора, возникающего в смешанной постановке и являющегося аналогом оператора Лапласа, который использовался в качестве переобуславливателя в итерационном методе сопряженных градиентов при решении эллиптического уравнения второго порядка с переменными коэффициентами. 4. Разработана методика распараллеливания двумерной и трехмерной задач для компьютерных систем, удовлетворяющих кластерной структуре, с использованием МРІ [65] - [68] технологии, которая может быть применена к другим задачам механики сплошной среды. Продемонстрированы результаты, показывающие высокую масштабируемость и эффективность алгоритма с точки зрения операций и обмена данными на многопроцессорных системах. 5. Создана вычислительная среда, включающая алгоритмы и программы обработки на многопроцессорных вычислительных системах сеточных данных, средства распределенного ввода-вывода. Данный программный комплекс позволяет при освоении нефтяного месторождения проводить анализ рисков и своевременно минимизировать их, совершенствовать технологии, обосновывать стратегические направления доразработки, улучшать показатели добычи, находить наилучшие интервалы вскрытия, исследовать поведение скважин и их групп, определять остаточные запасы, а также застойные зоны на конкретный период времени.
Архитектуры с распределенной и общей памятью
Из приведенных выше результатов, как и предполагалось, следует, что увеличение числа используемых ядер в какой-то момент приводит к замедлению роста производительности. Но даже в случае, когда число ядер максимально, то есть равняется минимальной размерности задачи, коэффициент ускорения больше 0.3р. В случае меньшего количества ядер программа демонстрирует хорошую масштабируемость. Коэффициент ускорения меняется в пределах 0.5р —;0.85р.
Изначально задача распараллеливалась именно под МРІ, т.е. под архитектуру с распределенной памятью. Однако с появлением многоядерных, многопроцессорных вычислительных узлов (в нашем случае это кластер нкс-зот, каждый узел которого это два процессора с четырьмя ядрами) архитектура кластера изменилась. Для передачи данных между узлами используется прямой доступ к памяти с помощью RDMA [69] технологии через сетевые карты InfiniBand [64], которые имеют низкую задержку (latency) и высокую пропускную способность (throughput). Однако память у процессов может быть общей, а значит пересылка данных между процессами, расположенными на одном узле будет происходить быстрее.
Все вычисления привиденные в работе производились с оптимальным заданием распределения процессов («многоядерньій»мрі). Например для четырех процессов выделялся один процессор с четырьмя ядрами, для восьми - один узел, состоящий из двух процессоров и четырех ядер. На кластере нкс-зот таким образом можно выделить до 1024 процессов, причем у каждых 8 процессов будет общая память. Если же рассмотреть «одноядерный» МРІ, ТО есть один процесс это один узел, то всего можно выделить только 128 процессов. На графиках рисунка 4.6 приведено сравнение коэффициентов ускорений на последовательностях сеток 64 х 8, 128 х 16, 256 х 32, 512 х 64, 1024 х 128 и 2048x256 с использованием до 64 узлов/ядер. Синий график ускорения, полученные для узлов («одноядерный» МРІ), а красный -ускорения, полученные для ядер («миогоядерньій»мрі).
Как видно из результатов, «одноядерный» МРІ более эффективен на небольшом числе процессов, когда вычислительная нагрузка достаточно велика и обмены данными не оказывают существенного влияния. На большом числе процессов распараллеливание с использованием ядер с общей памятью становится более эффективным.
В конечном счете, с ростом числа процессов коэффициент ускорения будет расти настолько быстро, насколько быстро пересылаются данные. Пересылка, данных на процессах с общей памятью происходит быстрее, поэтому коэффициенты ускорения «мпогоядерного» МРІ варианта растут быстрее, чем «одноядерного». Для нашей задачи данное отношение оказалось равным тройке. Ниже на рисунке 4.7 и в таблице 4.3 приводятся результаты, иллюстрирующие данный факт. 4,00
