Научная статья на тему 'Экспериментальное исследование некоторых разностных схем при сопряжении различных моделей фильтрации двухфазной жидкост'

Экспериментальное исследование некоторых разностных схем при сопряжении различных моделей фильтрации двухфазной жидкост Текст научной статьи по специальности «Математика»

CC BY
75
24
i Надоели баннеры? Вы всегда можете отключить рекламу.
Область наук

Аннотация научной статьи по математике, автор научной работы — Бочаров О. Б., Телегин И. Г.

Работа выполнена при частичной финансовой поддержке Российского фонда фундаментальных исследований, грант № 99-01-00622.

i Надоели баннеры? Вы всегда можете отключить рекламу.
iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

Experimental investigation of some difference schemes at conjugation of different models of two-phase fluid filtration

Numerical analysis of difference schemes for solving the problem on conjugation of different two-phase fluid filtration models is made. As a result of numerical experiments on different grids the optimal scheme is chosen.

Текст научной работы на тему «Экспериментальное исследование некоторых разностных схем при сопряжении различных моделей фильтрации двухфазной жидкост»

Вычислительные технологии

Том 7, № 1, 2002

ЭКСПЕРИМЕНТАЛЬНОЕ ИССЛЕДОВАНИЕ НЕКОТОРЫХ РАЗНОСТНЫХ СХЕМ ПРИ СОПРЯЖЕНИИ РАЗЛИЧНЫХ МОДЕЛЕЙ ФИЛЬТРАЦИИ ДВУХФАЗНОЙ ЖИДКОСТИ *

О.Б. БОЧАРОВ Институт водных и экологических проблем СО РАН Новосибирск, Россия e-mail: bob@ad-sbras.nsc.ru И. Г. ТЕЛЕГИН Горно-Алтайский государственный университет, Россия e-mail: dckt@bcentr.gorny.ru

Numerical analysis of difference schemes for solving the problem on conjugation of different two-phase fluid filtration models is made. As a result of numerical experiments on different grids the optimal scheme is chosen.

При изучении вытеснения нефти водой из пористых коллекторов с учетом капиллярных сил наиболее распространенной является модель Маскета — Леверетта (МЛ-модель). Однако, как показано А. Н. Коноваловым в [1], вследствие обращения в нуль функций, описывающих относительные фазовые проницаемости, естественные граничные условия для этой модели являются плохо обусловленными(градиенты решения в окрестности границы становятся бесконечно большими). В работе О. Б. Бочарова [2] предложено в качестве граничного условия на эксплуатационной скважине рассмотреть уравнение модели Ба-клея — Леверетта (БЛ-модель). В работах В. Н. Монахова, И. Г. Телегина [3, 4] эту модель предложено применять в окрестности эксплуатационной скважины. В итоге возникает задача сопряжения обеих моделей. Целью данной работы является сравнительный анализ некоторых численных методов решения данной задачи сопряжения.

Уравнения моделей. Одномерная МЛ-модель фильтрации двухфазной жидкости в однородной изотропной пористой среде в отсутствие массовых сил при заданном суммарном расходе фаз имеет вид [5]

msi = (k°a°(s)ai(s)sx — Q(t)b(s))w, (1)

где x £ [0,L] — пространственная координата; L — расстояние от нагнетательной скважины до эксплуатационной; t — время; s = (si — S°)/(1 — S° — S°) — динамическая насыщенность смачивающей фазы; si — истинная насыщенность смачивающей фазы;

* Работа выполнена при частичной финансовой поддержке Российского фонда фундаментальных исследований, грант №99-01-00622.

© О. Б. Бочаров, И. Г. Телегин, 2002.

) = const — остаточные водо- и нефтенасыщенности; m = m0(1 — S0 — S0), m0 — пористость коллектора; k0 = const — абсолютная проницаемость коллектора; a0(s) = fci k2, ai(s) = —pcs/(^2(k1 + ^k2)), pc(s) = (m0/fc0)1/27j(s) — капиллярное давление; 7 — коэффициент поверхностного натяжения; j(s) — функция Леверетта; b(s) = k1/(k1 + — коэффициент подвижности вытесняющей фазы; k = kj(s) — относительные фазовые проницаемости; ^ — вязкости фаз (индекс i =1 соответствует воде, i = 2 — нефти). Здесь и далее для функций fa = df/da.

При этом функциональные параметры модели имеют следующие свойства [5]:

1) 0 < (k1,k2) < 1 при sG(0,1), ao(0) = a0(l) = k1(0) = k2(1) = 0;

2) (fc1s,fc2s,pcs)gC(G), где G — замкнутая область в пространстве переменных (x, s), причем M—'Pcs)<M, M > 0;

3) | ln(a1),b/fc1|<Mo(M).

Полагая Q(t) = Q0 = const, введем безразмерные переменные: x = x/L, t = tQ0/(mL), при этом уравнение (1) запишется в виде

st = (ea(s)sx — b(s))x= — (2)

где e = y(m0fc0)1/2/(Q0L^2) — капиллярное число; a(s) = — a0js/(k1 + ^fc2); v = v(x,t) — скорость фильтрации вытесняющей фазы.

Значению параметра e = 0 соответствует модель Баклея — Леверетта. Обозначим через Г0 = {x,t|x = 0,t>0} нагнетательную скважину, через Гг = {x,t|x = М>0} — линию сопряжения моделей, через Г1 = {x,t|x = 1,t>0} — эксплуатационную скважину, а через Птга — область {x, t|m < x < n, t > 0}, m = 0; /, n = /; 1.

Сопряжение решений МЛ и БЛ-моделей с использованием смены эволюционной переменной. Рассмотрим случай, когда вблизи эксплуатационной скважины капиллярные силы малы. Разобьем область П01 на две, в области П0г фильтрация двухфазной жидкости описывается МЛ-моделью и соответственно этому в (2) e=0, а в области 0,ц — БЛ-моделью и в (2) e = 0. Предлагается для решения уравнения (2) в области П01 изучать следующую начально-краевую задачу:

s|x=0 = 1; s|t=0 = S0(x), x G [0,1]; [s]x= = [v]x= = 0, (3)

где [f]x= = {f(/ — 0,t) — f(/ + 0,t)} — скачок функции f(x,t) на линии Гг.

Из условий склейки на Гг следует условие easx|х=г-0 = 0, и задача (2), (3) распадается на две подзадачи. В области П0г

st = (easx — 6)л; s|t=0 = S0(x), s|x=0 = 1, easx|x= = 0. (4)

В области 0,ц фильтрация описывается БЛ-моделью, в соответствии с этим e = 0:

St = —bx; s|t=0 = S0(x), s|x=i = Sz(t). (5)

В работах [3, 4] предложено сменить эволюционную переменную в области 0,ц и искомую насыщенность s на расход v. В этом случае легче обосновать непрерывность расхода. Тогда задача (5) принимает вид

Vx = —^t(v); v|t=0 = v0(x), v|x=i = vi(t), (6)

где зависимость s = <^(v) является обратной по отношению к зависимости v = b(s), v0 = b[s0(x)], v(t) = b[S^(t)], функция S^(t) = s(x,t)|x= определяется по решению s(x,t) задачи (4).

Особенность задач (4), (6) заключается в склейке решений в = в(ж,£), (ж,£) € П0г уравнения (4), эволюционного по переменной ¿, и функции в(ж, ¿) = ^[^(х, ¿)], (ж, ¿) € Пц, выражающейся через решения ^(ж,^) уравнения (6), эволюционного по переменной ж.

Вопросам существования обобщенного решения задачи сопряжения МЛ- и БЛ-моделей (4), (6) посвящена работа [3], а первые предварительные расчеты проведены в [4]. Необходимо отметить, что в [3, 4] изучалась параболически регуляризованная задача для уравнения (6).

Построение алгоритмов численного решения задачи (4), (6) затруднено несколькими обстоятельствами:

1) заменой эволюционной переменной;

2) условием ^|^=0;1 = го, которое осложняет применение градиентных итерационных методов при линеаризации разностного уравнения.

Перейдем к описанию алгоритма решения задачи сопряжения (4), (6). Вначале произведем регуляризацию используемой в дальнейшем функции ^

( ^(в*), в€[0, в*],

^(в) = \ ^(в) [ ^(1 - в*), в€[1 - в*, 1].

Значение в* подбиралось по результатам численных экспериментов, черта над ^ (в) в дальнейшем изложении опускается.

Введем основную сетку в области П01 с распределенными узлами ш^т = (ж = гЛ, ¿п = пт; г = 0, У, п = 0,1, 2...} (Л = 1/У — шаг по пространственной координате, т = гЛ — шаг по временной переменной). Пусть при этом линии Гг соответствует номер Мг = /У узла на сетке (предполагается, что Мг не является дробным числом).

Также введем вспомогательную сетку в области Пг1 ш^ Тг, — (жг — гЛ^, tn = пт^; г = Мг, М1, п = 0,1, 2...} (Л^ = Г1Т^ — шаг по координате ж, т,и — шаг по времени), М1 соответствует линии Г1 или эксплуатационной скважине.

При этом г и г1 выбирались так, чтобы на линии Гг узлы сетки ш^^ совпали с узлами сетки ш^т.

Для нахождения численного решения использовались противопотоковые схемы из соображений необходимости вычислять в(ж,£) на линии склейки без привлечения дополнительных узлов сетки. В дальнейшем используются обозначения, принятые в [6].

Уравнение МЛ-модели в П0г аппроксимировалось с помощью неявной разностной схемы первого порядка

„п+1 _ „га р 'Ч___пп еп+1 „»

Л

г+ 1 г-1

^ _ ЬП++1, г =1,Мг _ 1, (7)

„п+1 _ „п 9 в м, в м 2е

Мг °Мг = _ ¿Ь „п+1 _ ьп+1 „п+1 = п =1

т = Л "мг -1 вх,мг ЬХ,Мг, в0 = в0 = 1-Для численного решения уравнения (5) использовалась разностная схема

„п+1 „п

I — в _

г L = 1 _ (1 _ ^п,г, г = Мг + 1,М1, вМ+1 = 57+1, (8)

т

где Ьп+1 линеаризовывались как Ьп+1 = Ьп + Ьпг(„п+1 _ „п). Схема (8) при а = 0.5 имеет погрешность аппроксимации О (Л + т2), при других значениях а € [0,1] погрешность аппроксимации О(Л + т).

Для численного решения систем (7), (8) применялся метод правой прогонки. Для аппроксимации уравнения (6) использовались два семейства разностных схем

_ „.п+1 ..га _ ..га ?.п+1 _ „.га _ „га а1-7--+ (1 _ -7- = _а2^г--(1 _ а2^г -' (9)

/¿V т,и

г = М ,М1 _ 1; <1 = 6(^Г+1);

1п+1 _ 1п+1 31п+1 _ + 7,п-1 31п+1 _ + 7,«-1

^г+1 1 _ ^ ,га+131г+1 41г+1 + ^г+1 п + 1 ПП^ -7- = _ а3^г -2--(1 _ -2-' ( )

» ¿V 2 / V 2 / V

г = М ,М1_ 1; 1М+1 = 6(^П+1)-

Схема (9) при а1 = а2 = 0.5 имеет погрешность аппроксимации 0(72 + т2), при а1 = 0.5 и 0"2 = 0.5 погрешность аппроксимации 0(7^ + т2), при а1 = 0.5 и а2 = 0.5 — 0(7^ + ), в остальных случаях погрешность аппроксимации 0(7^ + ).

Схема (10) при а3 = 0.5 аппроксимирует исходное уравнение со вторым порядком по обеим переменным, при а3 = 0.5 погрешность аппроксимации 0(7^ + т2). Схема (10) является трехслойной, поэтому на первом шаге, т. е. когда вода подходит к линии склейки, применялась двухслойная схема (9) с а1 = 1 и а2 = а3.

В численных расчетах использовался следующий набор параметров:

= в2; = (1 _з)2; 3 = (1 _в)/(0.9 + в); I = 0, 8; е = 0.5; ^ = 0.1; 5* = 0.05.

Выбор функциональных параметров был обусловлен простотой вычисления обратной функции ^(г>) = 6(-1)(^) для уравнения (6).

Ниже на рисунках жирными кривыми обозначены решения или характеристики, относящиеся к задачам сопряжения, а светлыми — результаты контрольного счета для БЛ-модели с а = 0.5 во всей области П01, пунктиром — линия Г.

Результаты расчетов. Первоначально тестировали разностные схемы на сопряжение смены эволюционной переменной, т. е. во всей области решали задачу для БЛ-модели без линии сопряжения и с ней.

Рассмотрим случай применения единой сетки ш^т на всей области. В численных расчетах использовались шаги сетки 7 = 7^ = 0.01 (Ж = 100), т = т^ = 0.001, М^ = 80, М1 = 100. Были проведены многовариантные расчеты при различных а, а1, а2, а3, з0(ж).

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

Сравнение численного решения задачи (4) - (6) с решением по БЛ-модели во всей области показывает, что после прихода воды на линию Г при малых з0(ж) имеет место дисбаланс порядка 2.3-2.4%, дисбаланс уменьшается с увеличением з0(ж) и /. Так, при 50(х) = 0.05 дисбаланс составил 0.6-0.7%. При з0(ж) = 0 фронт водонасыщенности по модели сопряжения движется быстрее. Для схемы (10) с а3 = 1 и I = 0.8 это приводит к уменьшению безводной нефтеотдачи на 6%. Еще одной особенностью решений задачи сопряжения является то, что при уменьшении а1, а2, а3 в схемах (9), (10) возникают осцилляции перед фронтом водонасыщенности. Они имеют место уже при а1 = 0.5, а2 = а3 = 0.72. С уменьшением з0(ж) осцилляции появляются даже при больших весах. Наилучшие результаты получились при решении задачи по схеме (10) при а3 = 0.74 и 50(х) = 0.05 (рис. 1).

Опыт использования разностных схем для задачи (4) - (6) на основной сетке показал, что оптимальной является более простая схема (9), так как, варьируя а1 и а2, можно

получить решение, близкое к тестовому. В то же время эта схема несколько хуже сохраняет баланс, нежели схема (10). Однако этот дисбаланс ложится в рамки допустимых пределов. График функции обводненности

1

„№ = ш>% /

0

для худшего варианта, когда дисбаланс достигает 2.45%, приведен на рис. 2. Интеграл в правой части уравнения вычислялся по формуле трапеций.

Теперь рассмотрим случай использования в области Пц вспомогательной сетки Тг,. В численных расчетах использовались шаги сетки Л = = 0.01, т = = 0.001 (Ж = 100), Мг = 80, Мх = 280.

Сравнение численного решения задачи (4) - (6) с решением по БЛ-модели во всей области показывает, что дисбаланс несколько меньше, чем при расчете задачи (4) - (6) на основной сетке, при в0(ж) =0 он составляет 1.8-2.6%, а при в0(ж) = 0.05 — порядка 0.29-1.0%.

0 0.2 0.4 0.6 0.8 х

Рис. 1.

Рис. 2.

О 0.2 0.4 0.6 0.8 X

Рис. 3.

Рис. 4.

При использовании схемы (9) с а < 0.57, а2 = 0.5 и в0(ж) = 0.05 имели место осцилляции после фронта водонасыщенности, при этом ухудшался баланс. С уменьшением в0(ж) осцилляции появляются даже при больших весах. Опять имело место опережающее движение фронта при в0(ж) = 0, как и на основной сетке. Расчет по схеме (10) при а3 = 1 и в0(ж) = 0.05 приводил к решению без осцилляций, но с нефизичной полочкой в окрестности фронтовой насыщенности. Наилучший результат получился при расчете по схеме (9)

при а1 = 0.58 и а2 = 0.5, з0(ж) = 0.05 (рис. 3).

Опыт использования вышеописанных разностных схем для задачи (4) - (6) на вспомогательной сетке показал, что оптимальной является схема (9).

На рис. 4 приведены результаты расчетов задачи (4)-(6) с е = 0.5, сравнение проводилось с решением по МЛ-модели во всей области П01. Расчет в области Пц проводился по схеме (9) при а1 = 0.58 и а2 = 0.5 на вспомогательной сетке. Вследствие условия еазх|х=г = 0 слева от Г образуется полочка. Капиллярные силы в МЛ-модели растягивают фронт водонасыщенности, поэтому время безводной нефтеотдачи по модели сопряжения больше, чем по МЛ-модели. Дисбаланс до прорыва воды составлял 0.6 %, а наибольший дисбаланс возникает после прорыва воды по модели сопряжения 1.75% при £ ^ 0.27. Такое расхождение вполне приемлемо, тем более, что с уменьшением е отличия в решениях уменьшаются, а в пластовых условиях обычно е < 0.1.

Список литературы

[1] Коновалов А. Н. Задачи фильтрации многофазной несжимаемой жидкости. Новосибирск: Наука. Сиб. отд-ние, 1988. 166 с.

[2] Бочаров О. Б. О задаче со сосредоточенной емкостью для одномерных уравнений двухфазной фильтраци // Динамика сплошной среды: Сб. науч. тр. / АН СССР. Сиб. отд-ние. Ин-т гидродинамики. 1985. Вып. 73. С. 149-155.

[3] МонАхов В. Н. Сопряжение основных математических моделей фильтрации двухфазной жидкости // Мат. моделирование. 2001 (в печати).

[4] Телегин И. Г. Численная реализация сопряжения основных моделей фильтрации двухфазной жидкости // Динамика сплошной среды: Сб. науч. тр. / РАН. Сиб. отд-ние. ИГиЛ. 2000. Вып. 116. С. 107-111.

[5] АнтонцЕв С. Н., КАжихов А. В., МонАхов В. Н. Краевые задачи механики неоднородных жидкостей. Новосибирск: Наука. Сиб. отд-ние, 1983. 316 с.

[6] Самарский А. А. Введение в теорию разностных схем. М.: Наука, 1971.

Поступила в редакцию 24 октября 2001 г.

i Надоели баннеры? Вы всегда можете отключить рекламу.