CC BY
42
УДК 519.6
А. С. Астракова 1, В. Н. Лапин 1, С. Г. Черный 1, О. П. Алексеенко 2
1 Институт вычислительных технологий СО РАН пр. Акад. Лаврентьева, 6, Новосибирск, 630090, Россия
2 Новосибирский технологический центр компании «Шлюмберже» ул. Зеленая Горка, 1/10, Новосибирск, 630090, Россия
Email: [email protected]
МОДЕЛЬ ФИЛЬТРАЦИИ ВЯЗКОПЛАСТИЧЕСКОЙ ЖИДКОСТИ В ЗАДАЧЕ ОПРЕДЕЛЕНИЯ ПАРАМЕТРОВ ТРЕЩИНОВАТО-ПОРИСТОЙ СРЕДЫ *
Построена модель фильтрации бурового раствора в трещиновато-пористую среду с вытеснением поровой жидкости. Для ее решения предложен оригинальный численный алгоритм, основанный на неявной конечно-разностной схеме. Обратная задача нахождения параметров трещиновато-пористой среды сформулирована в виде оптимизационной задачи. Реализованы два метода ее решения: метод золотого сечения и метод, базирующийся на генетическом алгоритме. Представлены результаты решения обратной задачи для различных групп варьируемых параметров.
Ключевые слова: модель радиальной фильтрации, вязкопластическая жидкость, трещиновато-пористая среда, вытеснение, обратная задача.
Введение
Процесс бурения скважины в породе может сопровождаться фильтрационными потерями бурового раствора. Одним из методов их устранения является технология закупоривания системы трещин в породе, находящейся в окрестности скважины, специально разработанными материалами (тампонаж). Успех применения такой технологии определяется правильным подбором тампонажного материала и режима его закачки и может быть повышен в случае, если известны параметры тампонируемой породы, такие как проницаемость, пористость, сжимаемость и др.
Целью настоящей работы является разработка методики определения фильтрационных параметров породы. Предполагается, что она имеет структуру специального вида - трещиновато-пористую. Трещиновато-пористая среда представляет собой совокупность пористых блоков, отделенных друг от друга трещинами [1]. Жидкость насыщает и проницаемые блоки и трещины. При этом размеры трещин значительно превосходят характерные размеры пор, так что проницаемость системы трещин кт значительно больше, чем проницаемость системы пор в блоках кП. В то же время трещины занимают гораздо меньший объем, чем поры, поэтому коэффициент трещиноватости тт - отношение объема, занятого трещинами, к общему объему
* Работа выполнена при финансовой поддержке Новосибирского технологического центра компании «Шлюмберже» и Интеграционного проекта СО РАН № 130.
Астракова А. С., Лапин В. Н., Черный С. Г., Алексеенко О. П. Модель фильтрации вязкопластической жидкости в задаче определения параметров трещиновато-пористой среды // Вестн. Новосиб. гос. ун-та. Серия: Информационные технологии. 2013. Т. 11, вып. 2. С. 18-35.
ISSN 1818-7900. Вестник НГУ. Серия: Информационные технологии. 2013. Том 11, выпуск 2 © А. С. Астракова, В. Н. Лапин, С. Г. Черный, О. П. Алексеенко, 2013
породы - существенно меньше пористости отдельных блоков тП. К трещиновато-пористым породам относятся некоторые виды известняков, песчаников, алевролитов, доломитов.
В основу предлагаемой методики положена модель плоскорадиальной фильтрации бурового раствора в трещиновато-пористую среду с вытеснением поровой жидкости. При этом буровой раствор представляется вязкопластической жидкостью, а поровая жидкость - ньютоновской. Уравнение пьезопроводности и закон Дарси, описывающие фильтрацию обеих жидкостей, имеют одинаковую структуру и отличаются только значениями входящих в них коэффициентов проницаемости и пористости. Поэтому процесс вытеснения буровым раствором поровой жидкости моделируется сквозным решением уравнений во всей области от скважины до удаленной на достаточное расстояние в область поровой жидкости границы с переключением значений коэффициентов уравнений на границе раздела жидкостей. Движение границы раздела жидкостей описывается отдельным уравнением. Заметим, что в постановке, когда буровой раствор и поровая жидкость предполагаются ньютоновскими жидкостями, задача была решена ранее (см., например, [2]). В случае фильтрации двух ньютоновских жидкостей значительно упрощаются вывод уравнений фильтрации и их численная реализация, а сама задача является частным случаем рассматриваемой в настоящей работе.
К конечно-разностной схеме решения уравнений пьезопроводности из-за разрывности их коэффициентов предъявляются особые требования, главными из которых являются обладание свойством консервативности и повышенным запасом устойчивости. В работе с целью создания численного алгоритма, обладающего указанными свойствами, строится оригинальная неявная конечно-разностная схема на основе дивергентной формы уравнения пьезопро-водности. Краевые условия для уравнений пьезопроводности - давление в скважине и в удаленной точке поровой жидкости. Решением этой задачи являются скорости фильтрации бурового раствора по трещиноватым и пористым блокам, по которым находится величина потерь бурового раствора.
Если из проводимых при бурении замеров известны временные зависимости давления в скважине и расхода потерь бурового раствора, то это позволяет поставить обратную задачу определения параметров трещиновато-пористой среды, сформулировав ее в виде оптимизационной задачи. В ней на наборе параметров трещиновато-пористой среды решается задача фильтрации (прямая задача) и минимизируется функционал отклонения между измеренной временной зависимостью потерь бурового раствора и зависимостью, предсказанной на основе разработанной модели фильтрации бурового раствора. Измеренная временная зависимость давления в скважине используется в качестве краевых условий для уравнений пьезопровод-ности модели фильтрации. Другими словами, суть методики заключается в подборе параметров трещиновато-пористой среды, которые обеспечат при решении прямой задачи фильтрации максимальное совпадение рассчитанной и замеренной временных зависимостей потерь бурового раствора. Подобранные параметры и являются решением обратной задачи фильтрации.
Модель радиальной фильтрации бурового раствора
в трещиновато-пористую среду с вытеснением поровой жидкости
Общие допущения. Фильтрация жидкости в трещиновато-пористой среде описывается в рамках модели плоскорадиальной фильтрации [1]. Жидкость фильтруется в радиальном направлении от скважины. Все параметры процесса зависят только от одной координаты - расстояния до скважины.
Предполагается, что проницаемости пористой среды и системы трещин по отдельности постоянны и равны кп и кт соответственно. Порода считается слабосжимаемой, т. е. пористости
породы тп и системы трещин тт линейно зависят от давления насыщающей жидкости р [1]:
тп = т0п (1 + Я (Р - Рпор )Х тт = т0т (1 + А (Р - Рпор где т0п и т0 т - пористости породы и системы трещин до начала бурения; вп и Д - коэффициенты сжимаемости пористой среды и системы трещин.
Буровой раствор полагается вязкопластической жидкостью, в которой касательное напряжение Т представляется в виде
Т = Т0 + К
с йм V2
йу
(1)
йм
где То - предел пластической деформации жидкости;--градиент скорости, перпендику-
йу
лярный направлению течения; К - коэффициент динамической вязкости, а п - степенной показатель в жидкости, подчиняющейся степенному закону.
Пористая среда пропитана поровой жидкостью. Поровая жидкость полагается ньютоновской жидкостью, в которой выражение для касательного напряжения Т является частным случаем выражения для напряжения в вязкопластической жидкости (1):
_ йм йу
Считается, что раствор вытесняет поровую жидкость, не смешиваясь с ней. Поэтому существует граница раздела жидкостей г = Яь (?).
Далее будут выписаны уравнения фильтрации вязкопластической жидкости в каждой из двух составляющих трещиновато-пористой среды, а также уравнение движения границы раздела жидкостей. Давление и скорость жидкости в пористой части обозначим через рп и мп, а в трещиноватой части - рт и . Из-за разницы давлений в трещинах и пористой части возникает переток жидкости из одной среды в другую. Величина перетока определяется коэффициентом перетока а0.
Уравнения пьезопроводности и законы Дарси модели фильтрации вязкопластической жидкости в трещиновато-пористую среду. Для вывода модели фильтрации вязкопластической жидкости в пористую среду рассмотрим установившееся плоскорадиальное течение вязкопла-стической жидкости между двумя параллельными пластинами, отстоящими на расстоянии Ж друг от друга. В этом случае уравнение количества движения превращается в следующее равенство [3]:
" / " ~лп ип С4п + 2V Тп
д-Р = -2 К С дг I п
Ж
п+1
п +1 ) Ж
Определим связь между действительной скоростью движения жидкости и и скоростью фильтрации м :
м = ти, (3)
где т - пористость, равная отношению объема пустот (пор) ко всему объему пористой среды. Пусть Н - высота слоя породы, в котором расположено N каналов ширины Ж. Тогда пористость данного пласта будет равна
NЖ
т-
Н
(4)
Из (2) имеем
и =
(2К)
1/п
Ж
1+1/ п
п
4п + 2
др +С4п + 2V Т дг
п +1 ) Ж
или, учитывая (3) и (4), получаем выражение для скорости фильтрации
м> = —-
1 К
(
N
пЖ
2+1/п
Л
1/п 1/п—1
Н (4п + 2)2"пК
др + С 4п + 2 V Т дг I п +1
_о_
Ж
1/п
которое будем трактовать как обобщенный закон Дарси:
м> = —
К
др дг
+ Т
1/п
(5)
В (7) введены обозначения для проницаемости k
k =
N
nW
2+1/n
1/П т^ 1/n-1
H (4n + 2)2"nK
и предельного касательного напряжения Т
Т =
4n + 2 Ï Т
(6)
(7)
п +1 ) Ж
Из (6) с учетом (4) выражаем ширину канала Ж через проницаемость к и пористость т :
1—2. _L f 4n + 2 k V+n W = K i+n 21+n f 4n + 2 '
n m
(8)
В свою очередь (8) используется в (7) для получения выражения Т от пористости т и проницаемости к пористой среды:
Т = Т0 — •0 K
1 1+n 11+n (4n + 2)1+nn1+n f m Vn
n +1
(9)
Таким образом, закон Дарси (5) с соотношением для Т (9) дают замкнутую систему для определения скорости фильтрации по давлению и параметрам пористой среды. Уравнение неразрывности в пористой среде [1] есть
д (рт) д (грлк)
- + -
= 0.
(10)
dt dr
Считая породу и вязкопластическую жидкость слабосжимаемыми,
M = шо(1 + ßnopoöu (Р - Рпор Ж
P = Po(1 + ßжидкости (Р - Ро)) ,
где р0 - плотность раствора при давлении p0р, ßp - коэффициент сжимаемости раствора, а
скорость фильтрации подчиняющейся закону Дарси (5), (9), из уравнения неразрывности (10) получим уравнение пьезопроводности для давления
dp 1 д
(
dt ß r dr
rk f dp
K \ dr
+ Т
1/n\
= 0,
где
ß ш0 ( ßпороды + ßжидкости ) .
Теперь мы можем выписать систему уравнений модели фильтрации вязкопластической жидкости в трещиновато-пористую среду.
Распространение давления в каждой из составляющих среды описывается своим уравнением пьезопропродности, при этом каждое из них учитывает фильтрацию как поровой жидкости, так и бурового раствора, отличаясь только значениями коэффициентов:
1 d
дРп___
dt ßj*r dr
f
rkn f dP>
1 d
dPj____
dt ßßr dr
K У dr rkт f dp
- + Т,
1/n\
+ -
an
K У dr
- + Т
J
1/n\
Kßß
- Рт )= 0,
a
Kßß
- Рт ) = 0,
(11) (12)
где
ßт = M0т (ßx + ßp ), ßn = ш0п (ßп +ß ) :
1—n
Т = Т
11+n 11+n (4n + 2)1+nn
nn
f m ^1+n
K
n +1
V ki J
, i = п, т.
n
Законы Дарси для скоростей фильтрации имеют вид
w = —
K
дРг_ dr
+ T
1/n
(13)
Входящее во все уравнения (11)—(13) выражение вида |--+Т \ вычисляется по сле-
дг
дующему правилу, отвечающему корректной передаче реологии неньютоновской жидкости,
dp
1/n
д-Р + т
dr
1/n
dp — + T dr
dp dr
— T
1/n
1/n
dp n — < 0 и dr
dp ~ — > 0 и dr
dp
dr dp
dr
>T,
>T
(17)
dp
dr
<T
Моделирование вытеснения поровой жидкости буровым раствором. Модель учитывает различия в реологических параметрах бурового раствора и поровой жидкости. Считается, что раствор вытесняет поровую жидкость, не смешиваясь с ней (рис. 1). Поэтому существует граница раздела сред г = Яь ) . В области между скважиной и границей находится буровой
раствор, за границей - поровая жидкость. Фильтрация обеих жидкостей описывается одними и теми же уравнениями. Учет среды производится путем переключения параметров
К, п, То в (11)—(13):
K, n,To =■
Кр, пр, T0 р > Rw < r < Rb,
к n t R < r
пор■> пор■> Опор' b ■
Рис. 1. Схема вытеснения поровой жидкости буровым раствором и краевые условия для уравнений пьезопроводности
Положение границы описывается задачей Коши
dR
dt
b _
= u(Rb,t), u = wT /тот, Rb|t=0 = R
well
Величина потерь бурового раствора из скважины Qwe^^ вычисляется по формуле
0ме11 = К + мт ) ,
где S = 2лRwell - площадь кольца единичной высоты на поверхности скважины.
Результаты определения параметров среды и реологии жидкости как для модели с вытеснением поровой жидкости буровым раствором, так и для модели без вытеснения, т. е. когда предполагается, что поровая жидкость имеет ту же реологию, что и буровой раствор, приводятся ниже.
0
Численный метод решения прямой задачи фильтрации бурового раствора в трещиновато-пористую среду с вытеснением поровой жидкости
Обобщенная запись уравнений пьезопроводности. Каждое из уравнений (11)—(13) может быть записано в следующем общем виде:
^ -дг дг
21 др + Т дг
+ Ьр - Ьр* = 0,
значения р, Б, 7, Т, Ь и р * приведены в табл. 1.
Значения коэффициентов уравнения (14)
(14)
Таблица 1
Уравнение р Б 7 Т Ь р *
(11) рп 1 вп*Г ГК К Тп ао кв рт
(12) рт 1 в*Г Гк Т К Т ао кв* рп
Неявная консервативная конечно-разностная схема для уравнения пьезопроводности. Уравнение (14) аппроксимируется неявной консервативной схемой с весами
т +1 т _ г-ч
р. - р. аО
Аг
{( - К+п ) + аЬ(р.+1 - р*т+1) -
(1 -а) Б 1.
{(+1/2 - ("/2) + (1 - а)Ь(р. - р]т) = 0, (15)
где
№т+1 = 7
"/+1/2 .+1/2
( рт+1 - рт+1
р+1 р +Т
А г Т+1/2 V .+1/2
\1/П] +1/2
/
"./+1/2
Г+1/2к К
, А = (А/+1/2Г + А-1/2Г)/2, А+1/2Г = Г/.+1 - Г/.
3+1/2
Верхний индекс т обозначает номер слоя по времени; нижний - / - номер узла по радиальной координате; а - параметр схемы. В схеме участвуют функции как из целочисленных узлов /, так и из дробных / +1/2, расположенных ровно посередине между узлами / и / +1.
При п ф 1 схема нелинейна и необходима ее линеаризация. Для этого на неявном т +1 слое вводится итерационный процесс, текущая итерация которого обозначается индексом 5:
{рт+1)
р5.
Итерационная схема записывается как
р5+1 - р5 аБ/
Аг А р] - рт (1 -а)Б
Аг
Представим
А
(5+1 = 7 /+1/2 /+1 / 2
Ш5+1 = 7 /-1/2 /-1/2
{(5+1/2 - (-1/2 ) + аЬ( р]+ - р*) = {(+1/2 -("1/2 ) + (1 -а) Ь( р. -р*т )
(16)
+ Т+1/2 )
(
( +1 +Т-1/2 )1/Л
где
+1 = р/+1 р/ +1 = р/ р/-1
А/+1/2Г
А/-1/2Г
Тогда линеаризация по Ньютону потоков Ж;!+1/2 даст
Ж-+1 = Ж- + 7
'' у+1/2 " у+1/2 ~ ^у +1/2
Ж*+1 = Ж* + 7
"у—1/2 ''у —1/2^^ у—1/2
п
— ( + Т../2 + )
у+1/2
п
—(
+Т—_„,) "н/'"' (+1 )
у—1/2
Введем обозначение и перепишем выражения для потоков (17):
$=р:+1
Ж-+1 = Ж- +_у+1/2
Жу+1/2 +1/2 + п
Ж-+1 = Ж- +
''у—1/2 ГГ]—1/2~
у+1/2^ у+1/2' 7
у—1/2 пу—1/2А у—1/2г
'Ту+1/2 (у+1 ),
Т—1/2 ( — —1 )
где
Т
] +1/2
= ( +Т+1/2 )1
1/п,—,,,— 1
Т —1/2 =(П +Т—1/2 ) ^ Подставляя формулы (18) в (16), получим следующий вид линеаризованной схемы:
П—1 — С А + С + р* + П + +1 = " у,
(17)
(18)
где
с7-В,7;
р = " у у+1/2 т С =_у у—1/2 Т 7
,0 1 У I 1 /О ч ^ -¿V 1 /1 ч ^
су7у—1/
п-+1/2А- Ау+1/2г
а
у+1/2' п п А- А гу—1/2' у+1/2 2
у —1/2 у—1/2 ^
— рт в
Р* ' ~Ау(ЖЖ+1/2 — Ж—1/2 ) + ЪРу — Ърр
гу .+. гу+1
V Ку Ку+1)
Аг А
Ж ±1/2 = Ж/2+(1—а)Ж/±1/2, =С+(1—а)рт, р* =ер;+(1—с)р*т.
Для правильной передачи реологии неньютоновской жидкости величину Т+1/2 в формулах (18) следует вычислять по правилу
Т
у+1/2
а)
б ) в)
Ау+1/2р < 0, А у +1/2р
А у +1/2г А у +1/2г
А у+1/2р > 0, А у+1/2 р
А у+1/2г А у+1/2г А у+1/2р
А у+1/2г
> Т+1/2 , то
>Т+1/2, то
А у +1/2 р +Т А г Т+1/2
у+1/2
А у +1^
пу+1/2
А у+1/2г
• —Т
у+1/2
< Т+1/2 , то 0,
а поток Ж]+1/2 в правой части О - по правилу
Ж,
]+1/2
а)
б )
б)
А ] +1/2 Р < 0, А ]+1/2 Р
А ]+1/2Г А ]+1/2Г
А ] +1/2 Р > 0, А ]+1/2Р
А ] +1/2Г А ]+1/2Г А]+1/2Р
А ]+1/2Г
>Т]+1/2, То +1/2
А ]+1/2 Р
А ] +1/2Г
-+т
]+1/2
>Т] +1/2, То ^]+1/2
<Т] +1/2 , То 0.
А ]+1/2 Р V А] +1/2Г
--Т
]+1/2
Построенная схема при <7 > 1/2 является абсолютно устойчивой и при ( = 1/2 аппроксимирует уравнение пьезопроводности с порядком аппроксимации О (Аt2 + Аг2), где
Аг = тах А ]+1/2г.
Скорость фильтрации рассчитывается из уравнения Дарси по известному распределению давления
т+1 1
w
]+1/2
']+1/2
-Ж,
т+1 ]+1/2 .
(19)
Решение уравнения для границы раздела бурового раствора и поровой жидкости Яь ^). Задача Коши (15) для границы раздела Яь ^) решается методом Эйлера 1-го порядка точности:
Ят+ = ят + Аt. и (Яьт, Г ), т = 0, 1,...
Я = ЯкеП, t0 = 0.
(20)
Поскольку скорость и в (20) определяется по скорости фильтрации бурового раствора в трещиноватой части породы, а та, в свою очередь, рассчитывается в дробных узлах (19), то
на каждом шаге метода (20) необходимо определять и (Яь, {п) в двигающейся точке г = Ят. Для этого сначала находится отрезок разбиения радиального направления, содержащий точку г = ят :
Ят е [Г г )
1УЬ с V]-1/2' ']+1/2/ •
Затем посредством линейной интерполяции по значениям скорости на концах этого отрезка вычисляется скорость в искомой точке г = ЯЯ^:
и (ЯЬ , Г ) = (ЯЬ - Г]-1/2)(и]+1/2 - и]-1/2)/( Г] +1/2 - Г]-1/2 ) + и]-1/2.
Верификация численного метода решения прямой задачи и анализ чувствительности решения
Плоскорадиальный фильтрационный поток упругой жидкости. Для верификации численного алгоритма использовалось аналитическое решение следующей задачи о неустановившейся плоскорадиальной фильтрации упругой жидкости, приведенное в [1]. В неограниченном горизонтальном пласте постоянной единичной толщины имеется добывающая скважина нулевого радиуса (точечный сток). Начальное пластовое давление во всем пласте одинаково и равно р^ . В момент t = 0 скважина пущена в эксплуатацию с постоянным расходом Qo . В пласте образуется неустановившаяся плоскорадиальная фильтрация упругой
жидкости. Распределение давления в пласте (в любой его точке и в любой момент времени) р(г, г) определяется уравнением пьезопроводности
* = Х А С г * \.
дt г дг V дг Начальные и граничные условия задачи следующие:
г = 0, 0 < г : р(г, г) = рк; г > 0, г = ~: р(г, г) = рк;
2пк С др
(21)
г > 0, г = 0:
Эг
= 00.
м v у г=0
Решением этой задачи является так называемая основная формула теории упругого режи-
ма фильтрации:
р(г, г) = Р(г, г) = рк
00М
4пк
-Е1
(
г
4хг
(22)
в которой введено обозначение для интегральной показательной функции
—Е1
( 2 ^ г
4хг
= I
е и
йи.
4хг
Для значений г2 / (4х-) < 1 с достаточной точностью выполняется
—Е1
(
г
4хг
Ы4^- — 0.5772
(23)
и в этом случае можно вместо (22) использовать формулу
р (г, г) = Р (г, г) = рк — °0М{ Ь^Хт — 05772\. (24)
4пкV г2 )
Из (22) может быть найден расход жидкости через цилиндрическую поверхность единичной высоты, радиусом г
0(г, г) = к др 2пг = 00е-г 2/4хг М дг
и скорость фильтрации
(25)
0,
2пг
0 „—г2/4х-
(26)
Поскольку для верификации алгоритма используется приближенная формула, оценим ее погрешность. Так, в [1] доказано, что при -е [10—'4,10 2] погрешность (23) не превосхо-
4хг
дит 5 %. Коэффициент пьезопроводности X в нашем случае имеет вид
к
Х = —^;. (27)
вК
Беря значения входящих в него параметров равными: к = 10—12м2, т0 = 0.3, впородЪ1 = 0, в
жидкости
4.58 10—10 Па"1, К = м = 10—3 Па - с,
получим X = 7.3 м / с.
В силу большого значения коэффициента пьезопроводности из (25), (27) следует, что стационарный процесс фильтрации достигается очень быстро на небольших расстояниях от
г
скважины. Кроме того, условие r2 / (4%t) < 1, при котором справедлива формула (24), выполняется в достаточно широких диапазонах изменения переменных t и r . Так, если время варьировать в пределах 60 c < t < 3000 c , а радиус - 1 м < r < 10 м, то значение величины
r2/(4Х) меняется в диапазоне [10-4,10-2]. Поэтому погрешность задания начальных и
граничных условий для (21) посредством формулы (24) составляет не больше 5% [1].
Для верификации предложенного в работе численного метода решения уравнения пьезопроводности он был применен к решению уравнения (21) при указанных выше значениях параметров в прямоугольнике [гл < r < rn ] X [tH < t < tK ]. На левой гл = 1 м и правой rn = 10 м границах задавались краевые условия, полученные из (24),
Р(Гл , t) = P(r,, t), p(Гп , t) = P(r„, t),
где P (r, t ) получено при значениях
pK = 37МПа, Q0 = 4.416 10—3м3 /c.
Начальное условие ставилось при tH = 60 c
Р(Г, tH ) = P(Г, tH X Гл < Г < rn.
Полученные при этих начальных и краевых условиях численные значения расхода и скорости фильтрации сравнивались с точными значениями, найденными из (25) и (26), вплоть до момента времени t = 3 000 c. Максимальная погрешность при этом не превышала 1 %.
Анализ чувствительности решения уравнений пьезопроводности. Если до этого мы решали линейное (n = 1 ) уравнение пьезопроводности для ньютоновской жидкости, то сейчас проиллюстрируем работу итерационной схемы и поведения решений уравнений пьезопро-водности в общем случае фильтрации вязкопластической жидкости в трещиновато-пористую среду с вытеснением поровой жидкости. Задавались следующие значения параметров модели
(11)—(13):
• для трещиновато-пористой среды
кп = 1 • 10—13 м2, Ш0п = 0.3, ßa = 0 Па
кт = 1.6 10—12 м2, т0т = 0.001, ßт = 0 Па—1;
• для вязкопластической жидкости (буровой раствор)
Кр = 3.83•Ю-2Па• с, np = 0.94, т0р = 4 Па/м, ßp = 4.58•1^10Яа"1^
• для ньютоновской жидкости (поровая жидкость)
Кпор = 10—2 Па • с, Пп0р = 1, Т0пор = 0 Па / м, ß^ = 4.58 •1^10Яа"1^
• коэффициент а перетока вязкопластической жидкости из одной упругой среды в другую и высота H слоя породы
а0 = 10—14, H = 10 м . Были выбраны следующие параметры схемы:
At = 0.5 c, min Ar, = 0.025 м, max Ar,- = 0.858 м.
J j
Начальные и граничные условия задачи:
t = 0, Рп(Rw, 0) = pwell(0), Рт (Rw, 0) = pwell(0), Rw < r < R^ : R^ = 0.1 м, R^ = 300м,
Рп(r, 0) = Рnоpовое, Рт (r, 0) = Рnоpовое, Рnоpовое = 37МПа
t > 0, r = Rw : Рп (Rw, t) = pwell(t), Рт (Rw, t) = pwell(t); t > a r = : Рп (R^, t ) = Рпоровое, Рт (, t ) = Рпоровое.
Измеренная зависимость pweu (t) представлена на рис. 2.
Рис. 2. Измеренная зависимость давления в скважине от времени
Характерные зависимости величин
^ = max
j
рП+ - p
П, j
^п, j
, = max
j
pi+j1 - p
т j
ps+l Гт, j
от итераций по s для первого и десятого шагов по времени приведены на рис. 3. Итерации для рис. 3 велись до тех пор, пока значения Err не достигали уровня точности представления чисел в компьютере. Во время решения практических задач точность сходимости итераций
для max } ограничивалась величиной 10 5.
Рис. 3. Кривые сходимости внутренних итераций на первом (кривые 1) и десятом (кривые 10) шагах по времени: сплошной линией отмечены профили ^ ; штриховой -
Для иллюстрации поведения решений обоих уравнений пьезопроводности полной модели фильтрации с вытеснением поровой жидкости на рис. 4 показаны распределения давлений в пористой (сплошная линия) и трещиноватой (штриховая линия) средах в различные моменты времени с нулевым значением коэффициента перетока (сверху)
a0 = 0
и ненулевым (снизу)
10
-14
a
0
39 -
ЗЭ ■
ЗЭ -
38 5 -
§
33 -
38 5 ■
§
33 ■
38.5 -
§
33 -
37.5 ■
37.5 ■
37 5 ■
37 -
10 20 30 40
г, m
37 ■
10
20 30 40
г, m
37 -,
10 20 30 40 г, m
39
38.5
§
33
37.5
37
10 20 30 40
г, m
39
38.5
43 §
33
37 5
37
10 20 30 40 г, m
Рис. 4. Радиальные распределения давлений рП (сплошная линия) и рТ (штриховая линия) в моменты времени г = 5 с (а), г = 500 с (б) и г = 4500 с (в) при а0 = 0 (сверху) и а0 = 10—14 (снизу)
Постановка обратной задачи определения параметров трещиновато-пористой среды и методы ее решения
Данные экспериментальных замеров. Имеются экспериментальные замеры давления
pwdl (t) в скважине и потерь бурового раствора QWfu (t) в зависимости от времени. Для
уменьшения шума на этих кривых проведено осреднение значений по пяти соседним точкам
1 2 _ 1 2 Qwell (ti ) = 5 ^ Qwell (ti+k ) > Pwell (ti ) = 5 ^ Pwell (ti + k ) ■
5 k=—2 5 k=—2
Зависимости до и после сглаживания показаны на рис. 5.
б
а
в
б
а
в
Рис. 5. Временные зависимости давления в скважине (слева) и потерь бурового раствора (справа), полученные из экспериментальных данных: тонкая линия - до сглаживания; толстая -
после сглаживания
F (x) =
(29)
Зависимость pwdl (t) используется для задания краевых условий для рп и рт
Рп (Rw,t) = Pwell(tX Рт (Rw,t) = Pwell(t) при решении уравнений пьезопроводности в прямой задаче фильтрации. Зависимость QwPu (t) является эталоном, к которому подгоняется расчетная зависимость потерь
QWm (t,x) путем подбора набора параметров трещиновато-пористой среды и вязкопластиче-ской жидкости
x = (x1,---,x16) = (кп ,т0п ,ß ,кт ,т0т ,ß ,KP, Пр, т0 р ,ßp ,Кпор, Ппор, Т0пор, ßnop ,a0,h) (28)
в модели фильтрации.
Оптимизационная постановка обратной задачи. Обратная задача определения параметров трещиновато-пористой среды формулируется в виде следующей оптимизационной задачи. Найти значения параметров x, обеспечивающие минимум функционалу
-T -|1/2
J (qZi (t) - о:ттг (t, x) )2 dt
_ 0 _
min F (x)
xeX
при фазовых ограничениях
X = {x: x'< x < <}.
Методы решения оптимизационной задачи. В ряде случаев бывает достаточно подобрать только один параметр среды, так как остальные слабо влияют на величину потерь бурового раствора. Как правило, таким параметром является проницаемость трещиноватой среды кт .
Поэтому в работе в таких случаях использовался метод решения оптимизационной задачи по одному параметру, а именно метод золотого сечения [4]. По сравнению с методами оптимизации по многим параметрам он обладает намного большим быстродействием.
Метод золотого сечения заключается в расчете значений функционала F(x) в нескольких точках на отрезке, выборе его наименьшего значения и сужении отрезка. Выбор точек производится таким образом, чтобы на следующем шаге можно было использовать значения функции в уже посчитанных точках. Метод состоит из следующих шагов.
1. Задаются границы интервала x0 = a , x3 = b и необходимая точность S.
2. Задаются точки на отрезке
Х3 Х0 Х3 Х0 , "^/З + 1
Х = Х3--—, Х2 = Х0 +-—, где Ф = —^ .
ф ф 2
3. Рассчитываются значения функционала в серединных точках
У = F (Х1), У 2 = F (Х2).
4. Из точек Х1, Х2 выбирается та, в которой значение функционала минимальное. На следующий шаг переносится она и ее соседи:
в случае yj > y2 берутся х0 = Х1, Х1 = Х2, y1 = y2 х2 = Х0 + -Х—— •
Ф
в случае у1 < у2 берутся х3 = х2, х2 = х, у2 = у1, х1 = х3 —0 .
ф
5. Рассчитывается значение функционала в изменившейся точке х1 или х2.
6. Шаги 4-5 повторяются до выполнения условия
--. < ю-з,
тах(| х1 |,| х21)
что обеспечивает поиск точки, в которой достигается минимум функционала с точностью до 0,1 %.
Можно заметить, что никогда не рассчитываются значения функции в границах интервала а, Ь . Но в этом нет необходимости, так как считается, что минимум функционала находится внутри отрезка.
Генетический алгоритм решения оптимизационной задачи по нескольким параметрам. В основу генетического алгоритма (ГА) положена теория биологической эволюции Дарвина, согласно которой популяция индивидуумов меняется в течение нескольких поколений путем селекции, рекомбинации и мутации, подчиняясь при этом законам естественного или искусственного отбора относительно предписанного критерия. Применительно к задаче определения параметров трещиновато-пористой среды при вытеснении поровой жидкости буровым растворам индивидуумом будет являться набор параметров х (28).
По сравнению с методом золотого сечения ГА позволяет эффективно находить решение задачи в случае варьирования нескольких параметров. Подробно об используемой версии алгоритма можно прочитать в монографии [5]. Схема алгоритма приведена на рис. 6 и характеризуется следующими шагами.
Рис. 6. Схема генетического алгоритма
1. Создание случайным образом начального приближения (поколения), состоящего из р различных конфигураций х0,...,х0,...,х0р - возможных решений задачи. Причем для каждого
индивидуума х0 = (х1,...,х16)0 выполняется х{ е [хЬ{,хК{] .
2. Вычисление значения функционала ¥ (х^) для каждого х£ . Определение лучших найденных решений хЬы ,1,..., хЬы р>^.
3. Если достигнута сходимость, то Выход, иначе переход на шаг 4.
4. Генерация нового поколения на основе хЬей1,..., х^ посредством операций ГА и переход на шаг 2.
Для построения новой популяции отбираются [Тг ■ р] лучших индивидуумов предыдущего поколения (этап селекции). Из них р раз случайным образом выбираются два индивидуума-родителя х = (х1,...,х16) и у = (у1,...,У16), по которым конструируется индивидуум
нового поколения те^ (и^,..., и16 ), где wi = х^ + Ц ( — xi), ai = гапй[—й, 1 + й] (этап рекомбинации). Каждый из р индивидуумов нового поколения слегка мутируется: иПеи = ±l(wLi — ^ )■ 2-16 гапй[0,1] (этап мутации). Здесь Тг, й, ¡1 - параметры селекции,
рекомбинации и мутации соответственно. Число поколений, необходимых для сходимости, обозначается ^^п. Значения параметров ГА, рассматриваемых в оптимизационных запусках, приведены ниже.
Результаты решения обратной задачи
Подбор параметра кт. Рассматривалась вариация параметра проницаемости системы трещин кт . Остальные параметры были зафиксированы на значениях, приведенных в анализе
чувствительности уравнений пьезопроводности. На рис. 7 представлена история поиска минимума функционала Я (29). Минимум функционала Я = 760,37 достигается при
кт = 1,61 -10-12 м2 / с и соответствует расходу 0^1 (?), приведенному на рис. 8 штриховой
линией.
Рис. 7. Поиск минимума Я методом золотого сечения
Рис. 8. Сравнение экспериментальной О^Хй (^) и рассчитанной ОСи (0 при оптимальном значении кт зависимостей
Подбор параметров кт, т0т, Кпор. Рассматривалась вариация параметров кт, т0т, Кпор в
диапазонах:
10-13м2 < кт < 3-10-12м2, 0.001 < т0т < 0.3, 0 Па • сп < Кпор < 0.0383 Па • сп
пор
Остальные тринадцать параметров были зафиксированы на значениях, приведенных в анализе чувствительности решения уравнений пьезопроводности.
Размер поколения Ы^еп составлял 200 индивидуумов, параметр селекции Тг = 0.1, параметр рекомбинации d = 0.7, параметр мутации ¡1 = 0.1.
История сходимости функционала Я (29) приведена на рис. 9.
Рис. 9. История сходимости функционала
Достигнутое минимальное значение функционала ¥ = 743.63 соответствует параметрам кт = 1.084 10—12 м2, т0т = 0.0129, Кпор = 0.00249 Па ■ сп. Зависимости значений каждого параметра на лучшем индивидууме поколения от номера поколения приведены на рис. 8. Сравнение экспериментальной О^р (Г) и расчетной ОИо? (Г) при оптимальных значениях кт, т0 т, К зависимостей представлено на рис. 10. Следует отметить, что при варьировании трех параметров кт,т0т,Кпор удалось достигнуть минимального значения ¥ = 743.63 , меньшего, чем ¥ = 760.37, полученного при варьировании одного параметра кт , что, однако, слабо сказалось на зависимости величины потерь от времени.
Рис. 10. Значения параметров кт, т0т, К на лучшем индивидууме поколения
Рис. 11. Сравнение экспериментальной 23 (Г)и рассчитанной ОИт (Г) при оптимальных значениях кт, т0т, Кпор зависимостей
Подбор параметров кт, т0п, а0, т0т, т0. Приводятся результаты решения обратной задачи с моделью фильтрации без вытеснения поровой жидкости. В этом случае полагается
Кр = Кпор = К, Пр = Ппор = П, Т0р = Т0пор = т0, вр = впор = Р.
Параметры варьировались в интервалах:
5-10—13м2 < кт < 3■Ю-12м2,
0.01 < т0т,т0п < 0.3,
—12
0 < а0 < 10 0 Па / м <т0 < 1 Па / м. Остальные параметры были зафиксированы на значениях:
кп = 1 10—13м2, вп = 0 Па—1, вт = 0 Па"1;
К = 3.83 10—2Па■ с, п = 0.94, р = 4.58 10—10Па_1; Н = 10 м.
—10
1.
Размер поколения еп составлял 300 индивидуумов, значения параметров ГА совпадают с приведенными в параметрах кт, т0т, Кпор. История сходимости функционала приведена на рис. 12. Оптимальное значение функционала Я = 764.28 соответствует параметрам к^ = 1.760 •Ю-12 м2, т0 г = 0.3, а0 = 0.787 •Ю-13, т0 = 0.01, Т0 = 0.185. Сравнение
экспериментальной О^еи (^) и расчетной ОСеа ) при оптимальных значениях кт, т0п, а0, т0т, Т0 зависимостей приведено на рис. 13.
0 10 20 30 40 0 2000 4000
с
Рис. 12. Ист°рия сходимости функци°нала Рис. 13. Сравнение экспериментальной ОЗ (0
и рассчитанной О^ (0 при оптимальных значениях кт, т0п, а0, т0т, т0 зависимостей
Для модели без вытеснения жидкости были также рассмотрены задачи с варьированием других групп параметров, решения которых приведены в табл. 2. Выделенные жирным шрифтом числа обозначают вариацию этих параметров в расчете.
Таблица 2
Оптимальные значения параметров и функционала
№ кр, м2 т0п а0 т0т т0, Па/м Я
1 1,106 -10 -12 0,3 0,200 -10 -12 0,01 0 851,71
2 1,106 -10 -12 0,3 0,105 -10 -12 0,01 0 850,67
3 1,013 -10 -12 0,3 0,175 -10 -12 0,3 0 845,21
4 1,812 -10 -12 0,3 0,737 -10 -13 0,01 0,2 764,80
5 1,580 -10 -12 0,3 0,200 -10 -12 0,256 0,047 766,98
6 1,760 -10 -12 0,3 0,787 -10 -13 0,01 0,185 764,28
Из табл. 2 видно, что при использовании модели без вытеснения не удается снизить значение функционала до значения, полученного по модели с вытеснением (функционалы принимают значения 764,28 и 743,63 соответственно). Более того, даже оптимизация по одному параметру по модели с вытеснением дает значения 760,37 меньшее, чем оптимизация по многим параметрам с моделью с единым описанием. Это позволяет сделать вывод о недостаточности модели с единым описанием бурового раствора и поровой жидкости для адекватного определения параметров среды. Судя по функционалам, модель с вытеснением поровой жидкости позволяет найти лучшее приближение к экспериментальному расходу и описать процесс фильтрации корректнее.
Заключение
В работе на основе уравнений неразрывности получены уравнения радиальной фильтрации вязкопластической жидкости в трещиновато-пористую среду, которые ранее в литературе не рассматривались. Построенная модель радиальной фильтрации в трещиновато-пористую среду также включает описание вытеснения поровой жидкости буровым раствором.
Для решения прямой задачи разработан эффективный численный алгоритм на основе неявной консервативной конечно-разностной схемы. При создании схемы учитывались требования на устойчивость и быстроту счета. При этом возникающие при построении схемы трудности были связаны со скачком значений коэффициентов сжимаемости, степенного показателя и предела пластической деформации на границе жидкостей, описанием ее движения, а также с необходимостью линеаризации уравнений.
Сформулирована обратная задача в виде оптимизационной и предложены надежные методы ее решения. Метод золотого сечения используется как наиболее быстрый, когда требуется подобрать лишь один параметр. В случае подбора нескольких параметров используется генетический алгоритм. Представлены результаты решения обратной задачи для различных групп варьируемых параметров. Исходя из оптимизационных расчетов, можно сказать, что модель с вытеснением поровой жидкости позволяет описать процесс фильтрации корректнее, чем модель без вытеснения поровой жидкости за счет нахождения лучшего приближения к зависимости измерений расхода от времени, полученной в результате эксперимента.
Список литературы
1. Басниев К. С., Кочина И. Н., Максимов В. М. Подземная гидромеханика. М.: Недра, 1993. 416 с.
2. Баренблатт Г. И., Ентов В. М., Рыжик В. М. Теория нестационарной фильтрации жидкости и газа. М.: Недра, 1972. 288 с.
3. Басниев К. С., Дмитриев Н. М., Розенберг Г. Д. Нефтегазовая гидромеханика: Учеб. пособие для вузов. М.; Ижевск: Ин-т компьютерных исследований, 2005. 544 с.
4. Каханер Д., Моулер К., Нэш С. Численные методы и математическое обеспечение: Пер. с англ. М.: Мир, 1998. 575 с.
5. Черный С. Г., Чирков Д. В., Лапин В. Н. и др. Численное моделирование течений в тур-бомашинах. Новосибирск: Наука, 2006. 202 с.
Материал поступил в редколлегию 04.12.2012
A. S. Astrakova, V. N. Lapin, S. G. Cherny, O. P. Alekseenko
MODEL OF VISCO-PLASTIC LIQUID FILTRATION IN PROBLEM OF FISSURED AND POROUS MEDIUM PARAMETERS DETERMINATION
In paper model of filtration of drilling agent into fissured and porous medium with extrusion of interstitial water was developed. To solve this problem original numerical algorithm based on implicit finite-difference scheme is proposed. Inverse problem of searching of fissured and porous medium's parameters is defined as optimization problem. Two methods of problem solution are implemented - method of golden section and method based on genetic algorithm. Results of inverse problem solving for different groups of vary parameters are presented.
Keywords: model of radial filtration, visco-plastic liquid, fissured and porous medium, extrusion, inverse problem.
