WWW.DISSERS.RU

БЕСПЛАТНАЯ ЭЛЕКТРОННАЯ БИБЛИОТЕКА

   Добро пожаловать!

Pages:     ||
|

f (x, y) dx dy = f (C), (1) hxhy C где (hx, hy) — размеры ячейки C. Эта функция может быть использована для консервативной аппроксимации уравнений решаемой на сетке задачи, содержащих пространственные производные до второго порядка включительно.

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

Разобьём ячейку на 4 интерполяционных прямоугольника P1, P2, P3, P4, построенных по девяти интерполяционным точкам p0,..., p8. Каждая точка pi (i = = 0,..., 8) содержит значения fi, dx, dy, которые аппроксимируют неизвестные знаi i n n n n 2 3 4 5 p p p 1 n n 1 6 P P 1 p p C 5 p n n 7 12 P P 3 p p 3 n n n n p 8 9 10 11 Рисунок 2: Малые соседи и интерполяционные точки чения f (pi), fx (pi) и fy (pi) соответственно. В пределах каждого прямоугольника величины fi, dx, dy интерполируются билинейным способом, формируя функi i (x, ции f y), fx (x, y) и fy (x, y).

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

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

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

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

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

i i Di = max fj - fk. (2) 0 j,k i В скалярном случае (когда fj — скалярные величины), вариация — это разность между максимальным и минимальным значениями интерполяционной функции в пределах ячейки:

( ) ( ) i i Di = max fj - min fj. (3) 0 j 8 0 j Цель алгоритма адаптации сетки — минимизация максимальной вариации среди листовых ячеек сетки при ограниченном числе этих ячеек:

max (Di) min при Nleafs Nneed, (4) iL где L — множество листовых ячеек, Nleafs — число листовых ячеек, Nneed — максимально разрешённое число листовых ячеек в сетке.

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

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

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

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

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

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

2. Если не требуется изменение числа ячеек в сетке (N = Nneed Nleafs), то алгоритм адаптации работает за время O (N) в среднем, и O (N ln (N)) в худшем случае.

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

Кроме того, возникает задача обобщения метода красно-чёрного упорядочивания узлов на адаптивную сетку.

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

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

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

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

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

c c c + vx + vy = 0, (5) t x y где c = c (x, y, t) — некоторая физическая величина, t — время, (vx, vy) = const — вектор скорости переноса величины c.

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

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

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

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

1 0.0.0.75 0.0.5 0.0.0.25 0.0 0 0.25 0.5 0.75 1 0 0.25 0.5 0.75 Рисунок 5: Результаты вычисления тестовой задачи на равномерной (слева) и адаптивной (справа) сетках. Сетки состоят из 4096 ячеек каждая.

Пунктиром показано начальное положение пятна Таблица 1: Среднеквадратичное отклонение результатов моделирования от точного решения 1024 ячейки 4096 ячеек 16384 ячейки Равномерная сетка 0,125 0,09 0,Адаптивная сетка 0,07 0,04 0,2.3 Вторая глава Во второй главе рассмотрена совместная фильтрация нефти и воды в пласте при нефтедобыче методом заводнения. Предполагается, что месторождение покрыто сетью водонагнетающих и продуктивных скважин. Через некоторые скважины в пласт вместе с водой поступает пассивная примесь — вещество, не влияющее на характеристики среды, в которой находится (например, индикатор).

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

.

.

.

.

.

.

.

.

.

.

.

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

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

q · F (S) — на источниках;

S m + div (F (S) K (S) gradP ) = Qw = (6) t q · F (S) — в остальной области;

div (K (S) gradP ) = q. (7) Здесь m — пористость; S — водонасыщенность; F (S) — функция Бакли-Леверетта;

K (S) — нелинейный коэффициент, включающий проницаемости и вязкости; P — давление; q — дебиты на скважинах; S — критическая водонасыщенность.

Скорость фильтрации флюида (закон Дарси): Wf = K (S) gradP ; скорость фильтрации воды: Ww = F (S) Wf. Нелинейные коэффициенты:

( ) kw (S) /µw kw (S) ko (S) F (S) =, K (S) = -k · +. (8) kw (S) /µw + ko (S) /µo µw µo Здесь kw (S) — относительная проницаемость воды; ko (S) — относительная проницаемость нефти; k — абсолютная проницаемость; µw – вязкость воды; µo — вязкость нефти.

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

(S · C + a) m + div (WwC + S (Ww)) = Qc = Qw · C. (9) t Здесь C — концентрация примеси в воде, a — концентрация примеси, адсорбированной на стенках пор (a = a (C) в случае обратимой адсорбции), Qc — источник 140 120 0.100 80 60 40 20 0 0 50 100 0 50 Рисунок 6: Рассчитанные поля водонасыщенности и концентрации примеси.

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

C Si = -Dij, i, j {x, y}. (10) xj Dij — эффективный тензор конвективной диффузии. Имеются феноменологические формулы для его получения из скорости Дарси и некоторых характеристик среды. В расчётах использована формула В.Н. Николаевского [4]:

Dij = [(1 - 2) ij + 2ninj] Ww, (11) где Ww = |Ww|, n = Ww/Ww, 1 и 2 — некоторые положительные коэффициенты. Можно видеть, что даже в изотропной среде имеется выделенное направление, определяемое вектором скорости Дарси.

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

.

.

.

.

.

.

.

.

.

.

.

.

.

.

.

.

.

.

.

= 6 · 10-13. Результаты моделирования (100 дней добычи) показаны на рисунке 6.

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

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

Pages:     ||
|



© 2011 www.dissers.ru - «Бесплатная электронная библиотека»

Материалы этого сайта размещены для ознакомления, все права принадлежат их авторам.
Если Вы не согласны с тем, что Ваш материал размещён на этом сайте, пожалуйста, напишите нам, мы в течении 1-2 рабочих дней удалим его.