CC BY
24УДК 622.276:531.72
МЕТОДИКА РАСЧЕТА ДВУХФАЗНОЙ ФИЛЬТРАЦИИ ЖИДКОСТИ В ПОРИСТЫХ СРЕДАХ НА ПРИЗМАТИЧЕСКИХ СЕТКАХ РАЗЛИЧНОГО ВИДА
О. И. Бутнев 1, В. Ю. Кузнецов х, В. А. Пронин х, М. Л. Сидоров х,
Р. М. Кац 2|, И. В. Афанаскин 2, А. В. Королев 2, П. В. Ялов 2
1 Федеральный научный центр Российский федеральный ядерный центр -Всероссийский научно-исследовательский институт экспериментальной физики,
г. Саров, virtualvi@mail.ru
2 Федеральный научный центр Научно-исследовательский институт системных исследований Российской академии наук», Ivan@Afanaskin.ru
В статье приводится методика расчета двухфазной фильтрации несмешивающихся жидкостей (нефти и воды) на основе модели «Black Oil» с использованием призматических структурированных и неструктурированных сеток. Для решения уравнений применяется полностью неявный конечно-объемный метод. Работоспособность методики проверена на задаче, основанной на стандартном тесте SPE8. Показано качественное и количественное согласие результатов с французским коммерческим симулятором Rubis компании Kappa Engineering.
Ключевые слова: двухфазная фильтрация, модель черной нефти, конечно-объемный метод, неструктурированные сетки.
TWO-PHASE FLUID FLOW CALCULATION METHODOLOGY FOR POROUS MEDIA
AT DIFFERENT PRISMATIC GRID TYPES
O. LButnevVV. Yu. Kuznetsov1, V. A. Pronin1, M. L. Sidorov1, R. M. Kats2|, I. V. Afanaskin2, A. V. Korolev2, P. V. Yalov2
1Russian Federal Nuclear Center, All-Russian Research Institute of Experimental Physics,
Sarov, virtualvi@mail.ru 2System Research Institute, Russian Academy of Sciences, Ivan@Afanaskin.ru
The article describes the two-phase fluid filtration calculation for non-mixing fluids based on the Black-oil model with an application of structured and unstructured prismatic grids. Equations solution is conducted through a fully implicit finite volume method. The described approach is tested on the SPE8 algorithm. Quantitative and quality parameters correlate with the commercial software package Rubis (Kappa Engineering).
Keywords: two-phase filtration, black-oil model, finite volume method, unstructured grids.
Введение. В настоящее время большая часть нефтяных месторождений на территории Российской Федерации находится на стадии истощения [1]. Доминирующей технологией разработки является заводнение (закачка в пласт больших объемов воды для вытеснения нефти и поддержания пластового давления). Поскольку в реальных пластах вытеснение нефти водой не является поршневым - для извлечения 30-35 % запасов нефти необходимо прокачать несколько поровых объемов воды. При этом в результате неоднородности пласта по фильтраци-онно-емкостным свойствам (по площади и по толщине), а также в результате существования обширной зоны совместной фильтрации нефти и воды, добывающие скважины длительное время одновременно добывают две жидкости - нефть и воду, причем доля воды в продукции скважин неуклонно растет. Для некоторых скважин соотношение дебитов воды и жидкости (обводненность) может составлять до 99 %. У многих скважин, находящихся в эксплуатации, обводненность составляет более 90 %. В связи со сложившейся ситуацией необходимо наличие программного обеспечения, позволяющего путем математического моделирования определять рациональные режимы разработки месторождений с целью увеличения коэффициента
нефтеотдачи и снижения производственных издержек.
На данный момент рынок программных продуктов, предназначенных для решения данных задач, представлен в основном иностранными компаниями: Schlumberger [2], CMG [3], Roxar [4], Kappa Engineering [5] и др. Однако доступ к данным программным продуктам в связи со сложной внешнеполитической обстановкой может быть существенно ограничен или полностью закрыт. Также при использовании иностранных программных продуктов, представляющих для пользователя «черный ящик», нельзя полностью исключить возможность утечки информации, представляющей коммерческую тайну.
Современное программное обеспечение для моделирования процессов разработки нефтяных месторождений должно удовлетворять ряду требований:
- наличие средств подготовки и анализа результатов (пре-пост-процессор);
- интуитивно понятный интерфейс;
- возможность проведения расчетов с использованием большого количества сеточных блоков (порядка 10 млн) и вычислительных ядер (до нескольких тысяч).
Геологическая структура нефтяных месторождений является достаточно сложной, имеющей неоднородности как по структурному составу, так и по геометрическим характеристикам (разломы, линзы и т. д.). Подробное воспроизведение в сеточной модели данных особенностей строения позволяет получать более точные результаты при прогнозировании процесса разработки и при оценке результатов геолого-технических мероприятий. Для описания подобных объектов необходимо построение специализированных сеток. Как правило, при моделировании нефтяных коллекторов используются сетки формата «Corner Point», которые являются регулярными топологическими и нерегулярными геометрическими. С ростом вычислительных ресурсов компьютеров стали все чаще использоваться неструктурированные сетки ввиду их способности хорошо описывать геометрические особенности модели (выклинивания, разломы и т. д.), поэтому представляются перспективными программная реализация и внедрение в практику повседневных расчетов использования таких сеток.
Модель двухфазной фильтрации
В большинстве случаев для моделирования разработки месторождений нефти на поздних этапах разработки можно использовать модель двухфазной фильтрации.
В данной работе мы опишем модель двухфазной фильтрации в системе вода - нефть, полученную на основе широко известной модели «Black Oil» [6].
Уравнения модели «нелетучей нефти» («Black Oil») представляют собой систему уравнений, объединяющих уравнения сохранения объема флюида в поверхностных условиях и уравнения движения Дарси. Движение среды описывается с точки зрения Эйлера. Система дифференциальных уравнений сохранения объема, описывающая двухфазное трехмерное течение жидкости в пористой среде, имеет вид:
dt
V Bo J
-V
Л
Wo
bo
VJ
= q0,
(1)
dt
V B* J
-V
Ww
Bw
= q*
где все обозначения - стандартные.
Уравнения движения (скорость фильтрации) описываются законом Дарси:
w0 = -Х 0 (Уро-у ^), (3)
^ = -Х,(Урт -у), (4) где Х0, Xкоэффициенты проводимости соответствующей фазы имеют вид:
кК , к К
Л. _ го а .Л _ гш а
х о =-, К =
Для замыкания системы уравнений (1)-(4) необходимы следующие соотношения:
Po - Pw = Pcow(Sw), (5)
So + Sw = 1. (6)
Отношение (5) выражает зависимость капиллярного давления в системе вода - нефть и используется для пересчета давления воды pw. Уравнение (6) используется для пересчета одной из насыщенностей.
Для полной формулировки математической модели исходная система уравнений (1)-(6) должна быть дополнена уравнениями состояния (данными PVT), а также начальными и граничными условиями.
Метод учета скважин в численных моделях фильтрации основан на допущении того, что вблизи скважины течение описывается аналитическим решением, граничные условия для которого определяются из численного решения задачи для пласта. Этот подход впервые предложен Писманом [7].
Скважина рассматривается как множество интервалов перфорации, вскрывающих сеточные блоки (ячейки). На скважине могут задаваться: забойное или устьевое давление, дебит флюидов в поверхностных или пластовых условиях. Распределение дебита по фазам и по интервалам перфорации осуществляется с учетом гравитационных сил и состава смеси в скважине.
В случае, когда на скважине заданы забойные давления, объемный фазовый дебит интервала перфорации, приведенный к стандартным условиям, имеет вид [8]:
qla= WI,\K.h S [p - pwb-y(Ht - Href )][; l = 1,2,.., Np ; * = o, w; (7)
l l [ ]J
qlg = q,oR,o; i = 1,2,..., np ; (8)
где p - давление в ячейке, содержащей интервал перфорации l; pw - забойное давление, приведенное к опорной глубине; WIt - коэффициент скважины для интервала перфорации l; h - толщина интервала l (может быть меньше, чем толщина ячейки, содержащей интервал l, но не может быть больше ее); St - индекс противотока; у - осредненный удельный вес жидкости в стволе скважины; Ht - глубина интервала l; Href - опорная глубина; N - число
интервалов перфорации.
Здесь для добывающей скважины должно выполняться соотношение:
Я =
1, Рг -Pwb -/(H -Href )> 0, 0, Рг - Pwb -У( H - Hf )< 0,
в котором для нагнетательной скважины знаки неравенств меняются на противоположные; у = ZZ qlayla/ZZ qia , причем ZZ qla * 0, в противном случае у = уа, где уа согласуется
l а / l a l а
с начальными условиями.
Коэффициент скважины (для прямоугольной сетки и сетки «Corner Point») в выражении (7) вычисляется как:
WI! =-
2ж
(\
ln
V rwl J
+S
где должно выполняться условие гы > тк1 ■ ехр(- ); гы - эквивалентный радиус сеточного
блока (ячейки), содержащего интервал 1; гк1 - радиус скважины на интервале 1; - скин-
фактор на интервале 1.
Эквивалентный радиус сеточного блока, содержащего интервал 1, для вертикальной скважины рассчитывается с помощью формулы Писмана:
гы = 0,28
4 i
k
yl
+ 4
V
yl
где кх1 и ку1 - абсолютные проницаемости сеточного блока, содержащего интервал 1, по соответствующим направлениям; Ах1 и Ау1 - размеры сеточного блока, содержащего интервал 1. Фазовые дебиты скважин рассчитываются с помощью выражения:
Ча = ^ Ч1а .
В случае, когда на скважине заданы фазовые дебиты ча, а = {о,н'} либо дебиты жидкости чь = ^ча , а ={о^}, используется упрощенный метод распределения дебитов по интер-
а
валам перфорации:
qla ~
BlaTla(la ; а = {o,w};
qla ~
Z BlaTl: l
BlaTlaqL
la
ZZ biPtIP
l p
; а = {o,w}; p = {o,w};
где Tla= WIk.hr fla
Bla№la
Kexp, аналогично для Tp .
lp •
Для случая, когда на скважине заданы фазовые дебиты либо дебиты жидкости, забойные давления рассчитываются из выражения (7).
Для случая, когда на скважине задан дебит жидкости, по формуле (9) производится распределение дебитов нефти и воды по интервалам перфорации. Затем по формуле (7), записанной для нефти или воды, определяется забойное давление. После чего по формуле (8) определяется дебит газа.
Искомые параметры скважины (забойное давление либо дебит) могут быть учтены в явной или неявной форме.
Под начальными условиями понимается начальное распределение искомых значений давления и насыщенностей на момент времени I = 0.
Под краевыми условиями понимается задание граничных условий (режимов работы) на границе моделируемой области и на каждом из источников, представляющих скважины.
В системе уравнений (1)-(6) 11 неизвестных:
Ро Р*
so Sw
Во в*
Ро Р*
к к ф
го Г* Т
для которых требуется 11 уравнений:
- 2 уравнения состояния или таблицы (PVT), позволяющие вычислить ро ,р* = f (Р) ;
- 2 уравнения или таблицы PVT, задающие Во, В* = f (Р) ;
- 2 таблицы относительных фазовых проницаемостей, задающих кго = f(Sw),
кг* = f (Sw );
- уравнение для насыщенностей, дающее So = 1 - Sw;
- таблица капиллярных давлений, дающая зависимость р* = ро - Рсо*),
Разностные уравнения полностью неявной схемы
Для получения разностных уравнений двухфазной модели фильтрации будет использован метод конечных объемов, а для аппроксимации потоков через боковые грани ячеек - метод отложенной коррекции.
Для решения разностных уравнений используется полностью неявный метод, известный еще и как метод совместного решения [6]. Суть данного метода заключается в записи производных от насыщенностей в правой части уравнений (1)-(2) через производные от давлений и решение уравнений в переменных давления нефти ро и насыщенности воды Sw.
Впервые это было предложено Дугласом [9].
Так как исходные уравнения являются нелинейными, в качестве метода линеаризации используется метод Ньютона - Рафсона [8]. В результате дискретизации и линеаризации образуются по два уравнения для каждого сеточного блока. При этом матрица принимает блочный вид.
Метод совместного решения. Обозначим через п временной шаг (п > 0). Введем следующий временной оператор для функции от времени V :
= vn+1 - V".
При неявной аппроксимации по времени система уравнений (1)-(2) определяется как
5 (_ Sy ^
at
v в; у
= А[^;+1(Ар;+1 -гГ'ай)] +
5 (фSc >
At
V Вс у
= А[ЛГ1(АрГ -^АЛ)] +
(10)
где At = tи+1 - Г . Система (10) является нелинейной с неизвестными р"а+1 и Sna+1, а ={с,;} и может быть линеаризована методом Ньютона - Рафсона. Для этого запишем:
и+и+1 „и+и , с,, гт+1,1+1 с*"+1,1 , со ~
Ра = Ра +5Ра' Sa = Sa а = С,; ■
где I обозначает номер итерации, а 5ра и 5Sa представляют изменение давления и насыщенности за этот итерационный шаг. Зная, что:
уИ+1 ^ у"+1,'+1 -."+1,'
= у"+и + 5у,
5у « у"+11 - у" + 5у
перепишем систему (10) в следующем виде:
1 А
/ \ "+1,1+1 V В; У
V В; У
+ 5
V В; У
= А[д;+и+1(Ар;+1,1+1 -гп;1М1т+%
а+1,1+1.
А
(_,А Г+1,1+1 кВс у
V Вс У
+5
V Вс у
= А^'+1 (Аро"+и+1 - Уа+и+1АИ)] + %
а+1,1+1.
(11)
В этой системе 5ра и 5Sa являются неизвестными. В качестве независимых переменных выбираются 5р, 55;.
Отметим, что 58а = -5Б;. Соответственно, левая часть системы (11) может быть представлена следующим образом: для водной фазы
5
V в; у
= с 5р + с „ 5S ;
где
с = ф_ сЯ
;р т Я
(S ^
V В; У
+
фSy
dв:
-Л'
dp
(ф1
V В; У
для нефтяной фазы
V вс у
= с 5р + с s 5S ;
Ср ± сЬс с'
где
Сор = фск
Во У
+
(1в-
-л1
dp
с я =
ГфЧ
Во У
Опуская далее временной индекс, запишем значение проводимости фазы на I + 1 ите-
рации
X =
I+1о I+1 Ма Ва
к,а = о, *,
где М+1 = М* . Таким образом, дебиты скважин определяются следующим образом:
Nw М*Т к'+1
€ -т) М((р$)м -РМ + Р'1 ) - *))*(* -х( 1т)
]=\ т=1
М
<£+1(',т) % ((р%>у+1 - рм-р+р %) - 7) х - х-т)).
]=1 т=1
Мо
Для удобства вводится понятие фазового потенциала
ф ? = р+1 - р'с+1 -у^+1vй, фо+1 = рм - у^Уй,
Представляя потенциалы, проводимости и дебиты фаз в основных переменных (8р и 8Sw ), получим.
Для водной фазы:
где
фм =Ф 1 + а 5р + й „ ss ,
* * * р г * Sw * '
V
d = 1 -
*р
V Ф У
dwSw =
^ V
dp
-Г (
V У
Для нефтяной фазы:
ф г+1 =ф1 + й 5р,
о о ор ± '
где
dop = 1 -
' dРo ^
V dp у
Проводимости фаз могут быть выражены аналогичным образом. Для водяной фазы:
Т*+1 = TW + Ewpëp + EwSwëSw,
где
Е =
Ак dB-\ ^
г*__
vМw Ф
К Е =
^к 1 ^
г*__
V М*в* у
К.
где
где
где
Для нефтяной фазы:
т:+1 = П + Е^др + Е^дБ*,
Еор =
к — го dp
( 1
кМоВо уу
К, ЕоБ* =
" ¿к, >
V
¿Б,
V V * у г^о о у
иВ
оо
К.
Дебиты скважин также представляются подобным образом. Для водяной фазы:
Ы* м*,
£1 = + ЕЕ™(',т) [е%др + + е|£дрй} ]д( х - х( -т))
'=1 т=1
е(= -
*р
1
е(') =
М М*
(
кг*
¿кг* ¿Б*
1 + ¿Р Р(¿Ы' - 7)
ар
*Рый
к1
кг*
М*
(РЫн) - р - Рсо* - р*р(¿Ын)))- К
¿Рс,
¿Б*
Для нефтяной фазы:
X* м
£1 = 1 + ЕЕ™(',т) \_еЦдр + еЩдБ* + еорындрЫЙ) ]д( х - х( >-т)),
'=1 т=1
,(') _
ор
= ^кг.
е(') =
¿р
(
Мо V
dМо ( р( ') - Р -РоР(НЫн) - н))-М
- ^ -р0Р(НЫН ) - н) )
1 + ¿Р М) - Н)
Ор
V
Мо V у
(')
Проделав подобные разложения для коэффициента объемного расширения Ва,а = о, *
получим
где
В1+1 = В1 (1 - Ы др), а = о,
а а ^ ар г /ъ ? ?
Ыар =-
' 1 ¿В а
V Ва ¿Р у
а = о,
Метод дискретизации уравнений на неструктурированных сетках
Аппроксимация по времени. Аппроксимация производной по времени искомой функции ф (немая переменная) в ячейке J по схеме первого порядка точности имеет вид:
/ \ п+1 п
[5ф| =фт -фт v ^ ^ тп
1
1
Аппроксимация по пространству. В случае, когда проводимость пласта задана в виде скалярной функции, значение градиента будем вычислять с использованием метода отложенной коррекции. В данном методе производная от давления по нормали в центре грани аппроксимируется с учетом вычисленного в ячейках сетки значения градиента, который определяется либо по формуле Гаусса, либо методом наименьших квадратов. При неявной аппроксимации значение градиента скалярного поля в ячейках определяется по значению поля с предыдущей итерацией (отложенной коррекцией). Отметим, что в случае использования разностных уравнений в дельта-форме аппроксимация потоков приращений через грани проводится без учета градиентов от приращений в центрах ячеек.
В задачах с анизотропной проводимостью градиент от искомой функции в центре грани определяется по формуле Гаусса интегрированием по контрольному объему, который является объединением среднемедиальных контрольных объемов относительно узлов граней ячеек.
Ранее подобный метод аппроксимации градиента был описан в работе [10]. Рассмотрим Р - ю ячейку и f - ю грань этой ячейки, которая является общей для Р и F ячеек. Точка f -центр грани. Проведем через центры Р и F ячеек плоскости, перпендикулярные нормали (иг)
к f - й грани. Найдем (рис. 1) точки Р' и F' пересечения этих плоскостей с нормальной линией к грани f .
Рис. 1. Шаблон точек для аппроксимации
v5n А
на внутренней грани
Векторы aRpp, и ARFF,, с учетом ориентации вектора нормали, определяются по следующим формулам:
Allpp = Rpi Rр =(Rf -Rр)- nf (nf • (Rf - Rp)), AR ee' = R E' - R E = (R f - R E)- n f (n f •( R f - R E)).
Считая, что в точках P' и F' известны значения ф , тогда поток через грань определяется по формуле:
(^grad Ф-n )f ASf =[a*b| ASf = (^)
9f' -Фр'
5n f f V Ф>f IAR
F'f
+
AR
p'f
ASf,
(12)
Значения фр, и фF, определим с использованием grad ф , вычисленного в центрах ячеек: Фр' = Фр +((grad ф)р •(Rp, - Rp )) = Фр +((grad ф)р ARpp), (13)
Фр' = ФР + ((grad Ф)Р • (- Яр )) = ФР + ((grad Ф)Р • АЯРР,). (14)
При подстановке (13)—(14) в (12), поток будет определяться по формуле:
^ 'А* =■
К),
АЯР'Г + АЯр,
(ф -Фр + ((ф)Р • А^рр' ) - ((gradф)Р • А^рр' ))
При неявной аппроксимации потока величины ф через грань , выделим неявную и явную части и запишем поток в виде:
v фэп у,
= хП (фР+1 -фР+1) + Я?,
где
Х? = 2
*(кгф)р +(1 -а)( кгф)р
КрКр
мф)р + кБф)Р ( Кр
АЯ
Р',
+ К
АЯ
р.,
А8,
Я? = x?
((§гаа ф)Р • аЯрр,) - ((§гадф)Р • аЯрр,)],
а, =
1 если (ро )р >(ро )р
0, если (Ро )р <( Ро )р
Пусть грань/принадлежит граничной поверхности, где задано значение функции ф :
Ф, = ФэG (,) ,
тогда поток через внешнюю грань будет определяться по формуле:
Эф
Ф, -фР'
АЯ
р,
АБ,.
(15)
Значение фр, вычисляется с использованием §гаё ф в центре ячейки (рис. 2) по формуле (13), при подстановке которой в (15), выражение для потока через внешнюю грань приобретает вид:
\фЭФ1 А* (фэо -Фр -((вгаёф)р -АЯрр,)).
V ФЭп У, 1 |АЯр,,
Выражение для потока через внешнюю грань при неявной аппроксимации, с учетом разделения на неявную и явную часть, здесь представим в виде:
v ф эп у,
а* =-(Хф), ФР" + Я?,
i)
где
Й =
(и
к
Мф)р
А8,
я! =х?
ФП^1 -((grad ф)Р-АЯ
рр'
Эф
Рис. 2. Шаблон точек для аппроксимации %- на внешней грани
Эп
Если грань / принадлежит поверхности, на которой задан поток, то
\ эфл
хфэП
V
ч ф( 1 ^)>
/ЯеЭО
тогда выражение для потока через внешнюю грань запишется в виде:
г
Эф
л Эпл
А^ =-Чф( 1 f,t
Для расчета значения градиента grad ф в ячейках разностной сетки используются два алгоритма. Первый алгоритм основан на вычислении значения grad ф в центре р - й ячейки по формуле Гаусса:
(grad ф)Р =
#фп^ Е фАА&1
^р__ 1=1_
аур
где ф1 - значение в центре грани (в точке 1 - й грани), п1 - внешняя нормаль к грани относительно р - й ячейки, - площадь грани, АУР - объем ячейки, суммирование ведется по всем граням ячейки.
Значение ф в центре внутренней грани определяется интерполяцией по длине кривой
)
фр дя к + фр дяи
дя и + дяк
Расчет тестовой задачи
Рассмотрим решения двухфазной задачи. Постановка данной задачи основывается на восьмом тесте из набора тестов SPE [11]. На рис. 3 показана начальная постановка теста SPE8. Проведем расчет данного теста на призматических сетках различного вида: прямоугольной (рис. 4) - сетка 1, многоугольной (рис. 5) - сетка 2 и треугольной (рис. 6) - сетка 3. Отметим, что сетки 2 и 3 имеют более разреженный вид в областях отсутствия внешних источников и стоков. Вдоль оси OZ использовались 4 слоя ячеек. В итоге сетка 1 имеет 400 ячеек, сетка 2 -108 ячеек, сетка 3 - 152 ячейки.
Рис. 3. Постановка теста SPE 8
Рис. 4. Ортогональная сетка (сетка 1)
Рис. 5. Многоугольная призматическая сетка (сетка 2)
Рис. 6. Треугольная призматическая (сетка 3)
В табл. 1-4 представлены свойства породы и жидкости для данной задачи. Область месторождения разрабатывается системой из трех скважин. Две добывающие скважины работают со значением дебита 2 000 м3/сут, нагнетательная скважина имеет приемистость 6 000 м3/сут.
Таблица 1
PVT свойства воды
Давление, атм Объемный коэффициент, м3/м3 Сжимаемость, 1/атм Вязкость, мПа-с
339,00 1,02 4,7E-05 0,36
Таблица 2 PVT свойства нефти
Давление, атм Объемный коэффициент, м3/м3 Расстворимость газа в нефти, м3/м3 Вязкость, мПа-с
339,00 1,55 0,0 0,397
Таблица 3
Талица относительной фазовой проницаемости для системы вода - нефть
SWAT KRW KROW PCOW
0,37 0 1 0
0,4046 0,013 0,242 0
0,4253 0,023 0,131 0
0,459 0,036 0,074 0
0,525 0,078 0,029 0
0,557 0,103 0,012 0
0,5839 0,143 0,01 0
0,6 0,191 0,008 0
0,63 0,502 0 0
1 1 0 0
Таблица 4
Таблица свойств коллектора
Слой Толщина, м Проницаемость, мД Пористость, д,ед,
1 7,62 500 0,3
2 22,86 50 0,2
3 22,86 20 0,2
4 45,72 10 0,1
На рис. 7-8 представлены графики дебитов добывающих скважин. Сравнение производится с коммерческим программным обеспечением (КПО) Rubis (Kappa Engineering) [5], расчет по которому осуществлен на ортогональной сетке. Из данных графиков видно, что время начала падения дебита нефти обеих добывающих скважин, рассчитанное по предлагаемой модели и с помощью КПО, хорошо согласуется друг с другом.
Рис. 7. Дебит нефти скважины «Ргойисег_1» Рис. 8. Дебит нефти скважины «Ргойисег_2»
На рис. 9 представлены временные зависимости значений насыщенности воды в ячейках со скважинами. Из графика видно, что данные величины имеют хорошее согласие результатов.
Рис. 9. Насыщенность воды в ячейках со скважинами
Далее проведем сравнение результатов, полученных на сетках 2 и 3. На рис. 10-11 представлены результаты дебитов добывающих скважин. Как видно из рис. 10, время начала падения дебита нефти на сетках 2 и 3 отличается от времени, полученного на сетке 1, примерно на 9 % для скважины «Producer_1». Для скважины «Producer_2» время начала падения дебита нефти для сеток 1 и 3 практически совпадает, в то время как для сетки 2 отличие составляет порядка 9 %.
2500
£ и
СП
^ 2000 О»
и 3 "С
£ 1500 Сч
3
В
в
ООО
я ей
а
и
В
N
•е-
3» 500
Н
«г ¡и
-Сетка 1 —Сетка 2 Сетка 3
200 400 600 800
Время (сутки)
1000
1200
1400
Рис. 10. Дебит нефти скважины «Рго(1исег_1»
2500
2000
£ и
Г-1
и1
«
и
■О
£ 1500
в
В
«
а
а о
я -
■ев и в
'С
и
ч:
1000
500
—Сетка 1 —Сетка 2 Сетка 3
200 400 600 800
Время (сутки)
1000
1200
1400
Рис. 11. Дебит нефти скважины «Producer_2»
На рис. 12 представлены временные зависимости насыщенности воды в ячейках со скважинами на разных сетках. Из рисунка видно, что на сетках 1 и 2 время прихода водного фронта к добывающим скважинам составляет одинаковую величину. Для сетки 3 отличие по времени составляет порядка 20 дней и возникает из-за различного количества ячеек между нагнетательной скважиной и добывающими скважинами «Producer_1» и «Producer_2».
На рис. 13-15 представлено поле насыщенности воды на момент окончания расчета.
Рис. 12. Насыщенность воды в ячейках со скважинами
Рис. 13. Поле насыщенности воды на момент окончания расчета (сетка 1)
Рис. 14. Поле насыщенности воды на момент окончания расчета (сетка 2)
Рис. 15. Поле насыщенности воды на момент окончания расчета (сетка 3)
Заключение. Создана методика расчета двухфазной фильтрации жидкостей на основе модели «Black Oil» с использованием призматических структурированных и неструктурированных сеток. Для решения уравнений применяется полностью неявный конечно-объемный метод.
Работоспособность методики проверена на различных тестах. В статье представлены результаты тестирования разработанной методики на тесте, основанном на стандартном тесте SPE8, и сравнение их с результатами расчетов с помощью зарубежного коммерческого симу-лятора Rubis (Kappa Engineering). На основе полученных результатов можно сделать вывод о работоспособности методики, так как достигается достаточное для практических расчетов качественное и количественное согласие результатов между разработанной методикой и Rubis.
Работа выполнена при поддержке гранта РФФИ 18-07-00671 А.
Литература
1. Методические указания по созданию постоянно действующих геолого-технологических моделей нефтяных и газонефтяных месторождений. Ч. 2. Фильтрационные модели. М. : ВНИИОЭНГ, 2003. 228 с.
2. Schlumberger : официал.сайт компании. URL: http://www.slb.ru/sis (дата обращения: 07.07.2018).
3. Computer Modelling Group Ltd. : официал. сайт компании. URL: http:// www.cmgroup.com/ software/stars.htm (дата обращения: 10.07.2018).
4. Roxar : официал. сайт компании. URL: http://www.roxar.com (дата обращения: 03.07.2018).
5. Kappa Engineering/Rubis - численное моделирование: описание. URL: https:// www.kappaeng.com/software/rubis/overview (дата обращения: 11.07.2018).
6. Aziz K., Odeh A. S. Comparison of Solution to a Three-Dimensional Black-Oil Reservoir Simulation Problem // JPT. 1981. January. Vol. 33. P. 13-25.
7. Peaceman D. W. Interpretation of well-block pressure in numerical reservoir simulation // SPE Journal. 1983. Vol. 23. № 3. P. 531-543.
8. Zgangxin C., Gyamrem H., Yuanle M. Computational methods for multiphase flows in porous media. Dallas : Texas : Southern Methodist University, 2006.
9. Douglas J. Jr., Peaceman D. W., Rachford H. H. Jr. A method for calculating multi-dimensional immiscible displacement // Trans. SPE AIME. 1959. № 216. P. 297-306.
10. Davidson L. A Pressure Correction Method for Unstructured Arbitrary Control Volumes // International Journal for Numerical Methods in Fluids. 1996. Vol. 22. № 4. P. 265-281.
11. Quandalle Ph. Eighth SPE comparative solution project: Gridding Techniques in Reservoir Simulation // SPE 25263. The 12th SPE Symposium on Reservoir Simulation. New Orleans : LA, 1993.
